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ABSTRACT 


The  equations  of  two-dimensional  motion  of  a  flexible  inextensible 
cable  are  linearized  by  a  small-perturbation  approximation  and  sinusoidal 
time  dependence  is  assumed.  The  simplified  equations  are  integrated  nu¬ 
merically  by  the  Kutta-Merson  method.  Separate  computer  programs, 

OMWAY  and  OMFLO,  have  been  written  for  the  Quadrant  1  and  Quadrant 
II  cable-towed-body  problems  and  are  listed  in  the  appendixes. 

ADMINISTRATIVE  INFORMATION 

This  research  was  supported  by  the  Naval  Air  Systems  Command  and  by  the 
General  Hydromechanics  Research  Program  under  SR  009  0101,  Task  0102  sponsored 
by  the  Naval  Ship  Systems  Command  Code  034 12B. 

INTRODUCTION 

Accurate  prediction  of  the  motion  of  cable-towed-body  systems  is  essential  for 
their  rational  design  and  for  simulation  of  their  operation.  For  the  simplest  of  these 
systems,  a  single  body  towed  by  a  cable  attached  to  a  towing  ship,  adequate  represen¬ 
tations  of  the  towpoint  motion,  the  dynamic  behavior  of  the  cable  as  it  is  affected 
by  the  tension  at  its  two  ends  and  by  hydrodynamic  and  gravitational  forces, 
and  the  motion  of  the  towed  body  are  needed.  If  only  steady-state  motion  is 
considered  where  the  entire  system  is  in  steady  rectilinear  translation,  analysis  is  possible 
by  numerical  integration  of  the  differential  equations  that  describe  the  cable  configur¬ 
ation  and  tension.  The  case  of  steady-state  motion  of  an  inextensible  tow  cable  in  the 
plane  described  by  gravity  and  the  direction  of  translation  has  been  solved  numerically 
by  Cuthill.1  The  tow  speed  and  the  tension  and  angle  at  the  towed-body  end  of  the 
cable  are  prescribed  and  the  viscous  force  on  the  cable  is  represented  by  functions  of 
cable  angle  based  on  experimental  data.  Integration  of  the  two  simultaneous  first-order 
nonlinear  ordinary  differential  equations  for  the  cable  tension  and  angle  as  functions  of 
distance  along  the  cable  is  performed  by  the  Kutta-Merson  method. 

The  unsteady  two-dimensional  motion  of  the  same  system  was  considered  by 
Whicker.3  Instead  of  prescribing  tension  and  angle  at  the  towed-body  end  of  the  cable, 
the  equations  of  motion  of  the  body  become  a  boundary  condition  for  the  equations 
of  motion  of  the  cable.  The  cable  equations  are  four  first-order  nonlinear  partial 


'Reference*  ire  lifted  on  pefe  46. 
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differential  equations  for  the  tension,  angle,  and  two  velocity  components  as  functions 
of  time  and  distance  along  the  cable.  Unless  the  restriction  to  an  inextensible  cable  is 
relaxed,  the  mathematical  classification  of  these  equations  is  parabolic  since  two  are 
hyperbolic  and  two  are  parabolic.  Consequently  they  cannot  be  solved  by  the  method 
of  characteristics. 

The  next  section  of  this  report  describes  the  system  of  inextensible  cable 
equations  given  in  Reference  2  and  a  particular  way  of  representing  the  viscous  force 
on  the  cable.  Then  two  important  simplifying  approximations  are  introduced:  the 
unsteady  dependent  variables  are  assumed  to  be  small  perturbations  of  the  steady-state 
variables,  and  their  time  dependence  is  assumed  to  be  sinusoidal.  The  first  of  these 
permits  the  linearization  of  the  unsteady  equations;  the  second  changes  them  from 
partial  to  ordinary  differential  equations  with  distance  along  the  cable  as  independent 
variable.  The  digital  computer  programs  used  to  integrate  the  equations  are  described 
in  the  next  section  and  are  listed  in  Appendixes  B  and  C. 

Thus  the  accuracy  of  the  analysis  reported  here  rests  on  the  validity  of  four 
main  assumptions:  (1)  that  the  motion  of  a  real  cable  can  be  represented  by  the 
solution  of  the  equations  of  motion  of  an  inextensible  cable,  (2)  that  it  remains  in  a 
vertical  plane,  (3)  that  it  consists  of  small  excursions  from  a  steady-state  configuration, 
and  (4)  that  these  excursions  are  sinusoidal  in  time. 


EQUATIONS  AND  BOUNDARY  CONDITIONS  FOR  TWO-DIMENSIONAL  CABLE  MOTION 

The  two-dimensional  motion  of  an  inextensible  cable  is  described  by  four 
nonlinear  first-order  partial  differential  equations  in  which  independent  variables  are 
distance  along  the  cable  s  and  time  /  and  dependent  variables  are  the  cable  angle  0, 
tension  T,  and  normal  and  tangential  velocity  components  U  and  V.  For  the  case 
where  the  motion  is  in  a  plane  containing  gravity  and  0  is  measured  from  the  hori¬ 
zontal  the  equations  are  those  derived  in  Appendix  A  which  are 
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where  ft  is  the  cable  mass  per  unit  length,  jytX  is  its  added  mass  (in  water)  per  unit 
length,  w  is  its  weight  in  water  per  unit  length,  and  Fi  and  Fq  are  the  normal  and 
tangential  components  of  the  viscous  force  per  unit  length  acting  on  the  cable.  These 
equations  were  also  derived  by  Whicker3  for  the  case  of  X=0. 

Two  cable-towed  systems  are  considered  in  this  report:  in  one  a  surface  ship 
in  a  seaway  is  towing  a  deeply  submerged  body,  and  in  the  other  a  deeply  submerged 
submarine  is  towing  a  float  which  is  slightly  beneath  the  free  surface  and  subject  to 
disturbances  from  its  seaway.  The  two  configurations  are  shown  in  Figures  1  and 
2.  In  either  configuration  the  (mean)  towing  direction  is  to  the  right  and  the  origin 
of  coordinates  is  at  the  center  of  mass  of  the  towed  body  with  x  in  the  (mean) 
towing  direction  and  y  up.  For  the  two  configurations  0  lies  between  0  and  x/2  and 
between  tr/2  and  ir  so  they  are  known  as  the  Quadrant  I  and  Quadrant  II  configurations 
respectively.  The  cable  is  of  length  L  and  the  quadrant  choice  of  $  causes  the  towpoint 
to  be  located  at  s~L  in  Quadrant  I  and  s=—L  in  Quadrant  II.  More  detailed  discus¬ 
sion  of  these  quadrant  conventions  is  provided  by  Springston.3 

In  the  Quadrant  I  problem  the  cable  motion  is  excited  by  the  effect  of  the 
seaway  on  the  towing  ship  so  the  boundary  conditions  at  s=L  are  the  motion  of  the 
towpoint;  the  boundary  conditions  at  the  bottom  are  provided  by  the  equations  of 
motion  of  the  body  tethered  by  the  cable.  In  the  Quadrant  II  problem  the  towpoint 
is  on  a  deeply  submerged  submarine  and  the  towed  body  is  subject  to  time-varying 
forces  due  to  the  seaway.  Therefore  the  boundary  conditions  at  s=0  are  the  equations 
of  motion  of  the  tethered  body  subject  to  these  forces;  the  other  boundary  conditions 
are  that  the  point  s=-~L  be  in  steady  rectilinear  translation.  Thus  each  problem  has 
kinematic  boundary  conditions  at  one  end  and  dynamic  conditions  at  the  other,  arid  in 
either  problem  the  excitation  is  due  to  the  action  of  the  seaway;  this  excitation  is 
through  the  kinematic  boundary  conditions  in  Quadrant  I  and  through  the  dynamic 
boundary  conditions  in  Quadrant  II. 

In  both  problems  the  steady  towing  speed  is  c.  In  the  Quadrant  I  problem  the 
towpoint  velocity  components  in  surge  and  heave,  positive  forward  and  up,  are  c  +  cs 
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(/)  and  Cf,  ( t )  and  are  assumed  to  be  known.  Therefore  the  kinematic  boundary  con¬ 
dition  for  Quadrant  1  is 

U(L,t)  sin  0(1,0  +  V{L,t)  cos  0  (1,0  *  c  +  cs  (/)  (5) 


and 


-  U(L,t)  cos  0  (L.t)  +  V(L,t)  sin  0  (L.t)  *  ch  ( t )  ^ 

In  the  Quadrant  II  problem  the  towpoint  is  in  steady  horizontal  translation  at  speed  c 
so  the  kinematic  boundary  condition  is 

U(-  L,t )  sin  0  (-  L.t)  +  V(-  L,t)  cos  0  (-  L,t)  =  c  (7) 


and 


-  C/(—  L.t)  cos  0  (-  L./)  +  K(- 1,/)  sin  0  (-  Z..Z)  =  0  <8> 

The  towed-body  equations  of  motion  used  here  are  those  of  a  very  simple 
representation  of  the  body  dynamics.  The  body  is  assumed  to  be  a  point  of  mass  m 
with  added  masses  in  surge  and  heave  m(  and  m*,  and  weight  (in  water)  WB,  on 
which  known  drag  and  lift  forces  Db  and  Lg  act.  No  body  pitching  is  represented, 
any  body  motion  which  is  excited  is  undamped  and  the  mass  center  of  the  body, 
center  of  buoyancy,  and  point  of  cable  attachment  are  coincident.  In  the  Quadrant  I 
problem  Dg  and  LB  are  taken  as  constants  but  in  the  Quadrant  II  problem  they 
represent  the  exciting  force  from  the  seaway  and  are  time-dependent.  If  UB  and 
VB  are  the  horizontal  and  vertical  components  of  the  body  velocity,  given  by 

UB(t)  ~  sin  0  (0,r)  +  V(0,t)  cos  0  (0,0 


and 


VB  (O  =  -  U(0,t)  cos  0  (0,0  +  K(0,0  sin  0  (0,0 
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differentiation  with  respect  to  /  and  substitution  into  the  body  equations  of  motion 
leads  to  the  dynamic  boundary  conditions 


I  dU  bV 

(m  +  ms )  <  -  (0,/)  sin  0  (0,/)  +  —  (0,/)  cos  0  (0,/) 

'  ot  ot  3^  \ 

+  [  £/(0,0  cos  0  (0,/)  -  K(0 ,0  sin  0  (0,/)]  (0,0  J  =  -£>a  ±  f  (0,/)  cos  0  (0,/)  (9) 


and 


(m  +mh )  (0,0  cos  0  (0,0  +  —  (0,/)  sin  0  (0,0 

+  [  U(0,t)  sin  0  (0,/)  +  K(0,0  cos0  (0,01  “(0,o|*  LB  -  WB  i  T  (0,0  sin  0  (0,0  (10) 

with  the  signs  in  front  of  the  T(o,t)  chosen  as  plus  for  the  Quadrant  I  problem  and 
minus  for  the  Quadrant  II  problem. 

In  steady-state  towing-cable  theory  the  viscous  force  components  Fj  and  Fq 
are  assumed  to  be  equal  to  D  times  known  functions  of  0,  where  D  is  the  drag  per 
unit  length  of  the  cable  when  0  is  x/2.  Fi/D  and  Fq/D  are  known  as  the  normal 
and  tangential  loading  functions.  In  this  report  it  is  assumed  that  the  viscous  force  on 
the  cable  in  unsteady  flow  can  be  represented  by  an  unsteady  generalization  of  the 
steady-state  force.  D  is  expressed  in  terms  of  Co,  the  coefficient  of  steady-state  drag 
per  unit  length  when  the  cable  is  normal  to  the  flow,  as  D=^y-  hCp  where  p  is  the 
water  density  and  h  is  the  cable  thickness,  or  dimension  in  the  direction  normal  to  the 
plane  in  which  the  motion  lies.  The  unsteady  generalization  of  this  expression  for  D 
consists  of  assuming  that  D  is  given  by 


D-  p 


U2  +  V2 
2 


hCD 


(11) 


where  Co  remains  the  drag  coefficient  for  the  cable  normal  to  a  steady  flow  but  U 
and  V  are  instantaneous  values  of  the  velocity  components.  The  unsteady  generali¬ 
zation  of  the  loading  functions  consists  of  replacing  the  functional  dependence  on  0 
with  the  same  functional  dependence  on  6,  defined  by  0=tan-1  U/V  since,  as  can  be 
seen  in  Figures  1  and  2,  this  is  the  angle  between  the  tangent  to  the  cable  and  the 
direction  of  cable  motion.  ( 8  clearly  reduces  to  0  when  U  and  V  take  on  their 

5 


J 


t 


l 


t 


I  I 


steady-state  values  of  U  -  csin  0  and  V  ±  c  cos  0.)  Comparison  of  Figure  2  with 
Figure  1  shows  that  the  positive  viscous  force  components  in  Quadrant  II,  Fr  and 

■  1  1  1  i 

-Fg  ,  should  exhibit  the  same  functional  depehdence  on  ir-0,  the  acute  angle  between 1 

i  f  ([  i 

the  cable  tangent  and  the  flow  direction,  that  Fl  and  Fq  exhibit  on  0  in  Quadrant  I. 

I  '  I  |  {  •  |  I  * 

Since  the  generalized  loading  functions  proposed  in  Reference  3  are  used  in  this  report, 
this  requirement  will  be  satisfied  if  the  loading  functions  are  given  by 


— -  =  A0^^  ±At^  cos0  +  cos  20  +  sin  6  ±B2  ^  sin  26 

D  1 1 


(12) 


and' 


Fc  * .  (nx/(n  ,  oiA  (n  ■  <n  ■  (d  . 

—  =  ±A0  +At  cos 0  ±42  cos  20  ±5,  sin0+fi2  sm  26 
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03) 


where  plus  signs  are  Used  for  Quadrant  I  and  minus  signs  for  Quadrant  II.  A  and  B 
are  constants  to  be  determined.  1 

i  1  .  i  ,  i 

An  example  of  how  A  and  B  may  be  related  to  experimental  data  is  provided 

I  •  *  1  *  i  | 

by  Eames,4  who  suggests  that  the  viscous  force  on  the  cable  can  be  assumed  to  consist 

i  i 

of  a  pressure  drag,  D(l-f)  sin2  0,  which  acts  normal  to  the  cable  tangent  and  a  , 
frictional  drag,  Df,  Which  acts  parallel  to  the  flow  direction  or  horizontal  for  the  steady- 
state  ( case.  Consequently  F £<  =  Z>(  1  -/)  sin2  0  +  Df  sin  0  and  Fq  t*  Df  cos  0,1  and  A 
and  B  are  given  by  Ao^  9  ,  A2^  -  -  ,  B,^  =  f,  *  f  and  the  1 

other  six  equal  to  zero.  This  representation  results  in  reasonable  values  at  the  ends'  of 

(  1  i 

the  range  of  0;  for  0  =  ir/2,:  Fl  is  D  and  Fq  is  zero,  for  0  =  0,  Fl  is  zero  and  Fg 
is  Df.  Furthermore  only  two  quantities,  Cp  and  •  /,  need  to  be  determined  from  ex¬ 
periments.  They  should  be  expected  to  be  dependent  on  Reynolds  number. 

The  parameter  X,  the  ratio  of  the  added  mass  per  unit  length  in  water  to  mass 

’  '  i 

per  unit  length  for  motion  in  the  x-y  pliane,  will  be  represented  in  the  form 


x-^» 

4/2 


(14) 
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Where  v  is  a  constant  to  be  determined.  Ideal-fluid  theory  predicts  *<*1  for  bare 
round  cables. 

In  either  quadrant  the  coordinates  x  and  y  of  a  point  on  the  cable  are  related 
to  5  and  0  by  -§£  =  cos  0  and  1  &  *  sin  0  so  they  are  given  by 
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I 


[15) 


and 

'  1  ,  ft 

>($,/)«  /  sin  0  (s'. 0  ds'  [16] 

,  ,  To  summarize,  U,V,  0  and  T  are  given  as  functions  of  s  and  t  by  solutions  of 

Equations  (1M4),  with  loading  functions  given  by  Equations  [12]  and  [13]  (with 
appropriate  signs),  D  given  by  Equation  [11],  X  given  by  Equation  [14],  kinematic 
boundary  conditions  provided  by  Equations  [5]  and  [6]  or  [7]  and  [8],  and  dynamic 
boundary'  conditions  by  [9]  and  [10]  (with  appropriate  signs).  Then  x  and  y  can  be 
found  using  Equations  [IS]  and  [16]. 

I.  I 

SMALL-PERTURBATION  FREQUENCY-DOMAIN 
'  EQUATIONS  AND  THEIR  SOLUTIONS 

Approximate  solutions  to  the  system  of  equations  given  in  the  previous  section 
are  obtained  for  the  case  when  the  exciting  disturbance  -  the  towpoint  motion  in  the 
Quadrant  I  problem  and  the  time-varying  force  in  the  Quadrant  II  problem  -  is  small 

»  f  s 

and  sinusoidal  in  time  by  assuming  that  0,  T,  U,  and  V  are  equal  to  the  steady-state 
values  plus  a  small  perturbation  term  which  is  proportional  to  the  cosine  of  hit  minus  a 
phase  angle,  to  is  constant  and  the  magnitude  and  phase  angle  of  the  perturbation 
terms  are  functions  of  s  found  as  solutions  of  eight  linear  ordinary  differential  equations 
obtained  from  the  four  nonlinear  partial  differential  equations,  [l]-[4].  Since  the 

J  l  !  '  * 

solutions  are  the  magnitudes  and  phases  of  a  sinusoidal  oscillation  they  are  known  as 
frequency-domain  solutions. 

The  steady-state  solutiqns  are  denoted  by  terms  with  subscript  zero.  For  ex- 

I 

ample  0(s,/)  is  given  by  1 

0U.O*  0oW  +  tor (J)  cos  [oit  -  60(5))  +  . . . 

=  0o(s)  +  0at  (r)  [cos  uit  cos  60(5)  +  sin  (tit  sin  60(5)]  +  . . .  U7] 

*  0o  (s)  +  0/j  (r)  cos  (tit  +  0/(j)  sin  (tit  + . . . 

=  0o (r)  +  fie  [0i(s)e-f«r]  + . . . 

.  I 
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where  Re  means  that  the  real  part  is  to  be  taken,  /  is  the  imaginary  unit,  and  +.  .  . 
indicates  that  0m  (s)  cos  [to/  -  6^(s)]  is  the  leading  term  in  a  perturbation  series  for 
0(s,/)  -  fa  C s).  <Pr  and  0/  are  related  to  the  magnitude  and  phase  0m  and  6^  by 


0*  *  0m  cos  6^ 


and 


0/  =  0M  Sin  60 


so 


and 


0m  2  +  0/2 


0* 


0,  is  given  by  0i  =  0/t  +  /0/. 

Similarly 

T\s.t)  =  Tots)  +  7m  (5)  cos  [to/  —  6r(s)]  +  . . 

=  70(s)  +  «e  [7, (*)<■•“*”]  +...  [18] 

with  corresponding  relations  between  7m,  &t,Tr,  7/,  and  7,.  As  was  seen  in  the 
previous  section,  C/0  and  V0  are  given  by 

U0(s)  =  c  sin  0O (s)  [19] 


and 


K0(s)  =c  cos  00  (5) 

and  satisfy  Equations  [1]  and  [2].  It  is  convenient  to  represent  U  in  the  form 


[20] 


U(s,t)  =  c  sin|0o(s)  +  0^^[0i ] }  +  fle[£/,(s)e'w']  +  . . . 

=  csin0o(5)  +  ccos0o(s)fle[0,(s)e'/a,/]  +  RelUxis)?10*1]  +  ... 


so  the  expansion  for  U  is 


LM 


U(s,t)  =  c  sin  0o(s)  +  fie  { [c0,  (s)  cos  0O($)  +  £/,  (s)Jeww' }  + 

In  the  same  way 

K(s,0  «  c  cos  0o(s)  +  de  { [—  c0,  (s)  sin  0O($)  +Vt  (s)]e'l(t3t  f  + . . . 
Then  if  Ux  -Ur  +  iUj  and  Uu  and  &u  are  given  by 


Um  -V(c^r  cos 0o  +  Ur)2  +  (C0/  cos 0o  +  Uj)2 


and 

C0/  cos  0o  + 

5(/  =  tan  1 - - - 

c4r  cos  0O  +  Ur 


Equation  [21]  reduces  to 

U(s,t)  -  U0(s)  +  UM (s)  cos  [ut  -  r  j(s)]  + . . . 

Similarly,  Vt  =  VR  +  iV/  and  VM  and  6k  are  given  by 


[21] 


[22] 


VM  =V(-c<Pr  sin  0O  +  Vr  )j  +  (-c0/  sin  0O  +  V/)2 


and 


*  A  -  c<t>i  sin  0o  +  Vj 

6k  *  tan  1 - ■ — : — - — 

-  c0*  sin  0o  +  VR 


to  give 

V(s,t )  -  V0(s)  +  VM(s)  cos  [wr-  6k(s)1  + .  •  • 

In  the  Quadrant  I  problem  the  excitation  is  due  to  speeds  in  surge  and  heave 
at  the  towpoint,  c,(t )  and  Cf,  (t),  yet  it  is  desirable  to  have  prescribed  displacements  in 
surge  and  hecve,  a,(t)  and  a*(/),  as  inputs. 

These  are  given  by 

at(J)sa,M  cos(wr-5,) 


4 


i; 


9 


and 


<ll,U)=UhM  cos  (CJ/-8*) 


so  i's  is  given  by 

<•,(/)  =  —  tofljw  sin  (cj/  —  8S) 

=  -  oxjSim  (sin  cj/  cos  6,  —  cos  cj/  sin  8,) 

=  csR  cos  co/  +  c,/  sin  c j/ 

=  Re  (c,ieltot) 

where  cs\  =  csR  +  ic,/  with  c,R  =  coa,*f  sin  8,  and  c,/  =  -cja,w  cos  8,. 
Similarly,  c/,  is  given  by 


[23] 


ch(t)s&e  (c* ieiU31) 

[24J 

where  c/,i  =  q,*  +  ichl,  chR  -  wa/,Af  sin  8h  ,  and  c*/  =  -0MhM  cos  8h 

The  excitation  for  the  Quadrant  II  problem  is  due  to  time-varying  drag  and  lift 
forces  on  the  towed  body.  Therefore  the  forces  are  given  by 


Dr  =  Dg  o  +  Du  cos  (cj/  —  8/> ) 
=  D,0+fle(/),e-,w') 


and 


Lb  ~  LBo  -  WR  +Lm  cos  (cj /  -  6*. ) 
=  Lbo  -  Kfc  +  (I,e/Wf) 


[251 


[261 


where  £>,  and  Z.t  are  zero  for  the  Quadrant  1  problem  and  given  by  Dx  =  £>*  +  iD/ 
and  I]  =  Lr  +  rX/  in  the  Quadrant  II  problem  where  DR  =  Dm  cos 
D/  -  Dm  sin  8q,  Lr  -  Lm  cos  8l  ,  and  Lj  -  Lm  sin  8l  . 

Substitution  of  Equations  [21]  and  [22]  into  Equation  [11]  gives 

D  =  —  { c2  +  2c  sin  0O  Re  [(c0,  cos  <p0  +  Ux  )c  't0f  1  +  . . . 

+  2c  cos  0o  K-  c'0i  sin  0O  +  V\  )<-,lut  1  +  . . . }  l-7l 
=  D0  { 1  +  --  Re  \cl(s)e^lu, ]  |  +  . . . 

where  Do  =  ( p/2)c2hCo  and  ct  =  U ,  sin  0o  +  cos  0o- 
Equations  [21]  and  [22]  give 
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jump 1  '"'wrran 


PCM 


U m  crin0o  +  fie  KC0,  cos0o  +  Ut)e‘lu>t]  ♦  . . . 
V  c  cos  0O  +  fle  l(—  c0 j  sin  0o  +  )«'/w/ 1  + . . . 


*  tan  0o 

1  +  fie  | 

V  c sm  0o  / 

1  +«e| 

r-c0,  sin  00  +  Vl  lait\k 

1  w  1  »  •  •  • 

=  tan  0o  | 

\ 

[l+«e| 

(  C  COS  0o  / 

/ c0j  cos  0o  +  Ut  ^  c0,  sin  0O  — 

1  L 

\  c  sin  0o  c  cos  0O 

Therefore 


6  ■  tan"'  i  tan  0O  +  tan  0O  fie  17—  C°S--+--1  +  —  **"**  —  ^)e  iui  1  + . . .  1 
I  \\  c  sin  0o  c  cos  0o  /  J  I 

and  use  of  the  Taylor  expansion  for  the  inverse  tangent  gives 

#.*„  +  ***•  «e  [~A»,  <=«♦.+ </.  +  *,  ifa».  -  K,  \  1  + 

1  +  tan3  0o  1_\  ^  sin  0o  c  cos  0O  /  J 

*0o+«e[0,(s)e-,w/j  +... 


[28) 


where 


$  m  tan  0o  /c0,  cos  0o  + 1/,  +  c0,  sin  0O  -  K,  \ 

1  1  +  tana0o  \  c  sin  0o  c  cos  0O  / 


U ,  K, 

=  0!  +  —  cos0o  -  — *■  sin  0o 
c  c 


Substitution  of  Equations  [27]  and  [28]  into  Equations  [12]  and  [13]  shows 
that  the  expansions  for  the  viscous-force  components  are 


Fl *  Fl  (Z>o,0o)  +2  Fl  (D0 ,0O )  fie 


dFL 
rf0  o 


(Do,0o)fle[0,(s)e-,wr]  +  ... 


[29] 


and 

FC  =Fg  (D0 ,0o ) +  2 Fg  (f?o,0o)«4-^e-,w'l  +-~^  (DOl0o)«e  [®,(»y'<*l +.. . 

L  c  J  o<p0  [30j 

The  two  steady-state  cable  equations  are  obtained  by  substituting  the  steady- 
state  terms  in  the  expansions  for  0,  T,  U,  V,  Fl,  and  FG  given  by  Equations  [17] 


MMMHK 
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through  (22],  [29],  and  (30]  into  Equations  [3]  and  [4]. 

Thus 

dd)n 

-  T0  — - Fl  (D0,4>o  )  +  w  cos  0O  =  0 

ds  [31] 


and 


dTo 

— - Fq  (D0.<Po  )  -  w  sin  0O  =  0 

as 


[32] 


Boundary  conditions  at  s  =  0  are  found  by  substitution  of  these  steady-state  terms  into 
Equations  [9]  and  [10].  The  result  is 

-DBo±T0  (0)  cos  0o  (0)  =  0 


and 

Lbo  —  Wb  ±  To  (0)  sin  0o(O)  =  0 


or 


To  (0)  =VDBo2+(LBO  ~  Wb) 


[33] 


and 


0o  (0)  =  tan"' 


LBo  ~~ 

—  Dbo 


[34] 


Dbo  is  positive  so  Equation  [34]  shows  that  Wb  must  be  greater  then  Lbo  for  the 
Quadrant  I  problem  and  less  than  Lbo  for  the  Quadrant  II  problem. 

Equations  [31]  and  [32],  which  are  simultaneous  nonlinear  first-order  ordinary 
differential  equations  for  0O  and  T0,  are  integrated  numerically  on  0  <  s  <  L  in 
Quadrant  I  and  0  >  s  >-L  in  Quadrant  II,  with  initial  conditions  provided  by 
Equations  [33]  and  [34].  The  integration  is  performed  by  the  same  Kutta-Merson 
subroutine,  KUTMER,  that  is  used  in  Reference  1.  Thus  the  results  of  Reference  1  are 
duplicated  although  it  is  necessary  to  write  separate  programs  for  Quadrants  I  and  II  and 
the  loading  functions  are  restricted  to  those  that  can  be  described  by  Equations  [12] 
and  [13]. 

Substitution  of  Equations  [17]  —  [22],  [29],  and  [30]  into  Equations  [1]  — 
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[4],  use  of  small-angle  approximations  to  the  multiple-angle  formulas,  subtraction  of 
Equations  [311  and  [32],  and  cancellation  of  a  mutual  factor  e'i0it  result  in 


dUx 

ds  1  ds 


/«0, 


[35] 


dtp  dVt 
‘  ds  ds 


»0 


[36] 


iun  [Xc0,  cos  0O  +  ( 1  +  X)  (/,  ] 


~To~~T'  ds  ~TFl  (D°’*o)c' 


and 


dFL 

d0o 


O>o,*o)0.  -  w0,  sin  0o 


[37] 


-i(A)H(Vt  +Xc0,  sin0o)«*~  -  —  Fc  (0o,0o)  c,  (0o,0o)  -  *0,  cos  0O  1381 

ds  c  a0  o 

Since  0O  and  T0  are  given  by  numerical  integration  of  Equations  [31]  and  [32]  and 
appear  as  coefficients  in  Equations  [35]  through  [38],  these  four  linear  ordinary  dif¬ 
ferential  equations  must  be  integrated  numerically.  This  is  also  done  by  the  subroutine 
KUTMER. 

Boundary  conditions  at  s  *  0  are  found  by  substitution  of  Equations  [17]— [22] 
into  [9]  and  [10].  Thus 

+  m,)[£/,(0)  sin0o  (0)+P,  (O)cos0o  (0)]=-Z),  ±  [7,(0)  cos0o  (0) 

-T0  (0)0,  (0)  sin  0O  (0)] 

and 

— /co(m  +m*)[  -£/i(O)cos0o  (0)+  K,  (O)sin0o  (0))  -LX±[TX  (O)sin0o  (0) 

+  Tp  (0)0,  (0)  cos  0o  (0)] 


or 

T,  (0) «  ±  (-  (0Xmr  -  mi, )  sin  0o(O)  cos0o(O)  +  K,  (0)[m  +  m,  cos20o  (0) 

+  mh  sin2  0o(O)] }  +  Dx  cos  0O  (0)  -  Lx  sin  0O  (0)) 
and 


[39] 
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NauriMI  mm*  i  *4 


140] 


0, (0)  =  *  *  ■  (- iu{-  Ux mm  +  mt  sinJ0o(O)  +  m*  cosa0o(O)]  -  Vx(0 )(m,  -  w* ) 

*o  (0) 

sin0o(O)cos0o(O)}  —  Dx  sin0o(O)  —  Lx  cos0o(O)) 

with  plus  signs  for  Quadrant  I  and  minus  signs  for  Quadrant  II.  Recall  that  Dx  and 
Lx  are  zero  for  the  Quadrant  I  problem. 

Boundary  conditions  at  s  =  L  for  the  Quadrant  1  problem  are  obtained  by 
substitution  of  Equations  [17]  and  [21]  -  [24]  into  Equations  [5]  and  [6].  Thus 

Ux(L)-cs\  sin  0o  (L)  —  Ch\  cos 0o  (L)  [41] 


and 


VX(L)  ■  c,i  cos  0o  (D  +  Chi  sin^o  (^)  [42] 

Substitution  of  Equations  [17],  [21],  and  [22]  into  Equations  [7]  and  [8]  shows  that 
the  boundary  conditions  at  s  -  -L  in  the  Quadrant  II  problem  are 


£/,(-/.)  =  0 


[43] 


and 


Vx(—  L)  =  0  [44] 

Satisfaction  of  the  two-point  boundary  conditions  is  achieved  by  assigning  differ¬ 
ent  sets  of  initial  values  to  the  dependent  variables  at  the  lower,  or  passive,  end  of  the 
cable  and  integrating  Equations  [35]  -[38]  up  the  cable  twice,  once  with  each  set  of 
initial  values,  so  that  two  different  solutions  are  obtained  along  the  cable.  Then  the 
linearity  and  homogeneity  of  the  lower  boundary  conditions  and  the  equations  permit 
linear  superposition  of  the  two  solutions  such  that  the  inhomogeneous  boundary  condi¬ 
tions  at  the  upper,  or  excited,  end  are  satisfied.  This  is  known  as  a  “shooting”  method 
of  solving  a  two-point  boundary-value  problem. 

The  two  sets  of  solutions  are  called  the  A-mode  and  B-mode  solutions  and  the 
solution  to  the  full  problem  is  given  by 
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=  CA  UxA  +Cb  U\b 


[45] 


I 


Vi‘CA  VxA  +CB  VlB  [46] 

0i  =Q0i/<  +Q0ifl  [47] 

and 

Tx*CATlA  +CBTlB  [48] 

where  the  subscripts  on  Ux,  Vx,  0,,  and  Tx  denote  the  solution  mode  and  CA  and  CB 
are  constants  to  be  determined  below. 

In  the  Quadrant  I  problem,  values  of  UxA (0),  UlB( 0),  VxA (0)  and  K,fl(0)  are 
assigned  and  TlA  (0),  TxB(0),  $lA  (0)  and  tlB(0)  are  found  from  Equations  [39]  and 
[40]  with  plus  signs  and  Dx  •  Lx  =  0.  Next  Equations  [35] -[38]  are  integrated  on 
0  <  s  <  L  so  that  both  modal  solutions  are  known  on  0  <  j  <  I.  Then  Equations 
[45],  [46],  [41],  and  [42]  show  that  CA  and  CB  must  be  given  values  such  that 

CAUlA  (L)  +  CBUlB(L)=  c,i  sin0o(£.)  -  c*i  cos^o(f-) 
and 

CA  V,A  (L)  +  CB  VxB(L )  =  c,i  co^o(i)  +  Chi  sin0o(f>) 
are  satisfied.  Therefore 

c  _  [c»i  sin0o(£)-c<n  cos0o(f,)l  VxB(L)-[ctX  cos <t>o(L)+  chi  sin0o(Z.)l  UiB(L)  49j 
A  UlA(L)VxBa)-U.B(L)VxA(L) 


1 

1 

J 


and 


_  [c,i  cos0p(f»)  +  Cfi i  sin0o(I)]  UxA(L)-[  ct\  sin0o(I)-c»i  cos0o(£)]  VxA(L)  [5Qj 
B~  UxA(L)VxB(L)-UxB(L)VxA(L) 


i 
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Finally  Ux,  Vx,  0,,  and  Tx  are  found  on  0  <  i  <  L  by  substituting  back  into 
Equations  (45) -[48]. 

The  Quadrant  11  problem  is  solved  in  much  the  same  way.  Values  are  assigned 
to  <I>ia  (mL),  0i /»(-£),  TlA(-L),  and  TxB(-L).  Then  Equations  [43] and  [44]  show  that  UlA (-L), 
UlB(-L),  VxA(-L)  and  ViB(-L)  are  all  zero.  (Otherwise  Equations  [45]  and  [46],  with 
left-hand  sides  equal  to  zero,  show  that  CA  and  CB  would  be  determined  by  conditions 
at  s  *  -L).  Equadtions  [35] -[38]  an  integrated  on  -Z,  <  s  <  0  so  that  both  modal 
solutions  are  known  on  -L  <  s  <  0.  Equations  [45] -[48]  and  [39]  and  [40]  with 
minus  signs  show  that,  if  P  and  Q  are  defined  by 

P  =  T,  (0)  —  icj  {  Ux  (0)(ms  —  mh )  cos0o(O)  sin0o(O)  +  K,  (0)[m  +  ms  t;os20o(O) 

+  mh  sin20o(O)]  ( 


aiiu  I 

Q=  To(O)0i(O)-iw  |f/i(0)[m  +  m  ,  sin20o(O)  +  m/,  cos20o(O)] 

-  F,(0)(m,  -  mh)  cos0o(O)  sin0o(O)| 

with  PA,  PB,  QA,  and  QB  defined  with  appropriate  modal  values  of  Tlt  0,,  Ux,  and 
Vx  on  the  right-hand  side,  the  boundary'  conditions  at  s  =  0  become 

CAPA  +  CBPB  S—Dx  cos0o(O)  +  £,  sin0o(O) 


and 

CaQa  +CbQb=Dx  sin0o(O )  +  /.,  cos0o(O) 

Therefore  CA  and  CB  are  given  by 

_  [-Dx  cos0o(O)  +  Lx  sin0o(O)]  QB-[DX  sin0o(O)  +  Ii  cos0o(O)] />g  [51] 

A  PaQb-PbQa 


and 


_  \D\  sin0o(O)  +  Lx  cos0o(O)]  PA  -  [-Dx  cos<p0(0)  +  Lx  sin0o(O)]  QA  [52] 
B.  P4Qb-PbQa 

Finally  Ux,  K, ,  0, ,  and  Tx  are  found  on  -L  <  s  <  0  by  substituting  back  into  Equation 
[45]-[48]. 
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The  analysis  is  completed  by  computing  the  cable  configuration.  This  is  denoted 
by  x  and  y  and  is  found  from  Equations  [15]  and  [16].  Equation  [17]  suggests  that 
x  and  y  are  given  by  expansions  of  the  form 

x(s,t)  *  x0(s)  +  fie  [xiCs^'^']  +  . . .  [53] 


and 


y  ( s.t )  =  y0(s)  +  fie  ly ,  (s)e,,w'l  +  . . . 


[54] 


Then  if  xt  *  x*  +  ix/,  xm  =  Vxr2  +  r,-1 ,  and  6X  =  tan"1 -  it  is  found  that 

*r 

x  *  x0  +  xm  cos  (w/  —  6X)  +  .  .  . 

Similarly  yx  =  yR  +  iyi,  yM  +  yi2,  and  8y  =  tan  1  — 

yR 

Substitution  of  Equations  [17],  [53],  and  [54]  into  Equations  [15]  and  [16] 
and  use  of  the  small-angle  approximations  to  the  multiple-angle  formulas  give 


*o 


i: 


cos0o  (s')  ds' 


[55] 


X,  =  -P  0,(s')  sindo (s')ds' 


[56] 


y  o  *  s:  sin0o(s’)cfo' 


[57] 


and 


yi  *  f:  0,(s')  cos^o  (s')dj' 


[58] 


These  integrals  are  to  be  evaluated  for  0  <  s  <  L  in  the  Quadrant  I  problem  and  for 
0  >  s  >  ~L  in  the  Quadrant  II  problem.  Because  of  the  choice  of  origin,  the  con¬ 
figuration  given  is  that  seen  by  an  observer  moving  with  the  towed  body. 

After  the  complex-valued  functions  0,.  7*,,  c0,  cos  0O  +  U\>  <0i  sin  <f>0  + 

Vx,  X|,  and  yx  have  been  found  the  amplitudes  and  phases  of  the  perturbation 
quantities  are  easily  found. 

There  is  one  redundancy  in  this  small-perturbation  frequency-domain  solution  to 
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the  two-dimensional  cable  problem  since  there  is  no  way  of  assigning  a  meaning  to  the 

i  •  t  , 

time  t  =  0.  Consequently  the  phase  angle  of  either  the  horizontal  or  vertical  ex- 
citation  is  meaningless  and  may  be  set  equal  to  zero. 


COMPUTED  PROGRAMS 

ii  ,  / 

The  digital  computer  programs  which  solve  the  Quadrant  I  and  Quadrant  II 
problems  are  OMWAY  and  OMFLO  respectively,  and  are  listed  in  ( Appendixes  B  and  C, 

*  1  i 

together  with  sample  calculations.  They  are  written  in  FORTRAN  IV  for  a  CDC  6700 
computer.  ' 

>  Although  Equations  [35]— [38],  (S3],  [541 ,  and  the  boundary  conditions  are 
written  for  complex-valued  functidns  of  s,  they  must  be  separated  into  real  and 
imaginary  parts  because  the  integrating  subroutine  KUTMERi  is  not  available  in  complex 
arithmetic.  However,  the  existence  of  a  complex -conjugate  relationship  between  the 
real  arid  imaginary  pahs  makes  solution  for  two  instead  of  four  modes  sufficient. 

Whenever  possible  the  notation  of  the  previous  section  is  retained;  some  important  ex- 

I  1  :  !  : 

ceptions  are  the  angle  0  which  is  replaced  by  PH,  and  the  frequency  u>  which  is 

1  1  I  , 

replaced  by  OM.  In  general,  quantities  which  pertain  to  the  tqwed  body  start  With  a  B. . 

!  ,•  i  1 

Thus  m, ,  the  body  added  mass  in  surge,  is  BMS.  The  cable  is  of  thickness  CH  and  is 
N-\  segments  long;  each  segment  is  of  length  DS.  Its  massi  per  unit  length  n  is 

s'  I  '  / 

replaced  by  ULM  and  its  added  mass  parameter  v  by  AMP.  Its  weight  per  unit  length  , 
in  water  is  WUL  and  the  density  of  water  is  DN.  Further  explanation  of  the  program 
is  provided  by  comment  cards  in  the  listing.  <  , 

i  Each  program  includes  four  subroutines.  AT  A  prevents  the  inverse-tangent 

■  I  /  ^  i  j 

function  from  dividing  by  a  prohibitively  small  number  if  an  angle  is  very  close  to  > 

±  it  12.  KUTMER  performs  integrations  with  a  fourth-order  Kutta-Merson  method.  This 

i  ’  i  >  1  '  ,  f  1 

is  a  generalization  of  the  Runge-Kutta  method  that  provides  an  automatic  reduction  in 

! 

integration-step  size  when  an  error  criterion  is  qot  met.  DAUX  and  DARN  provide 

integrands  for  KUTMER;  DAUX  for  0O,  T0,  x0,  and  y0  and  DARN  for  Ux,  Vt,  0,, 

1 

and  Tt.  KUTMER  is  not  used  for  x,  and  y ,;  they  are  calculated  by  a  simple  first-  < 
order  integration  within  the  mam  prograni.  , 

»  1  ,  1  i  l  .  I 

In  the  program  OMWAY,  J  is  an  integer  which  increases  from  /  =  1  at  the  towed 
body  (s=0)  to '  J-N  at  the  towpoint  {s=L).  Sqccessive  data  cards  are  used  to  specify; 

•  i  1  . 

1.  Name  of  program,  OMWAY.  (  i 

2.  Body  lift  BL,  drag  BD,  and  weight  BW. 

3.  Magnitude  and  phase  of  surge  and  heave  of  towpoint,  ASM,  DAS,  AHM,  DAH, 


i  I 


i 


t 


and  frequency  OM. 

4.  Body  mass  and  added  masses  in  surge  and  heave,  BM,  BMS,  and  BMH. 

5.  WUL,  ULM,  DS,  and  N. 

1  6.  Tow  speed  C,  DN,  CH,  cable  drag  coefficient  CD,  and  AMP. 

,  7.  Constants  for  Equation  [12]  to  give  Fl/D0. 

8.  Constants  for  Equation  [13]  to  give  Fc/D0. 

!  , 

These  data  are  written  out.  Do  is  computed  and  Equation  [14]  is  used  to  compute  X, 

which  is  represented  by  AMC/ULM.  Equations  [33]  and  [34]  are  used  to  find  To 

. 

and  fa  at  J=l  and  Uo,  Vo,  <t>o>  To,  xo,  and  yo  are  computed  by  using  Equations  [12], 

[13],  [19],'  [20], i  [31],  [32],  [55],  and  [57]  and  are  written  out  as  functions  of  J. 

Values  are  assigned  at  7*1  as  UjA  real  and  V]A  zero,  <t>iA  and  TjA  are  found  from 
Equations  [39]  and  [40],  and  the  A-mode  solution  is  computed  by  integrating 
Equations  [35] -[38]  on  1  <  J  <  N.  Then  Uxb  is  set  zero  and  VlB  is  set  real  at 
/*  1  land  the  5-mode  solution  is  found  in  the  same  way.  Next  ctX  and  chi  are  com¬ 
puted  and  the  required  values  of  Ux  and  Vx  at  J*N  are  found  from  Equations  [41] 
and  [42]  and  CA  and  Q  ate  found  from  Equations  [49]  and  [50].  Ux,  Vx ,  0,  and 
r,  are  found  along  the  cattle  from  Equations  [45] -[48]  and  their  magnitudes  and 
phases  are  computed  and  written  out  Finally  Equations  [56]  and  [58]  are  used  to 
find  xx  and  yx  and  their  magnitudes  and  phases  are  computed  and  written  out. 

OMFLO  uses  the  nomenclature  of  OMWAY  wherever  possible.  J  is  an  integer 
which  decreases  from  J=-l  at  the  towed  body  (j  =  0)  to  J=-N  at  the  towpoint  (s=-L). 
For  performing  integrations  down  the  cable,  which  is  done  when  <f>o.  To,  xo,  yo,  *i, 
and  yj  are  computed,1  an  integer  K  is  used  which  is  defined  by  K=-J  and  thus  in¬ 
creases  from  K=1  at  the  towed  body  to  K-N  at  the  towpoint.  Since  these  inte¬ 
grations  are  made  in  the  direction  of  decreasing  s,  the  integrands  computed  in  DAUX 
are  negatives  of  those  in  QMWAY.  U],  Vj,  4>j ,  and  7/  are  computed  by  integration 
in  the  direction  of  increasing  s  so  an  integer  I  is  defined  by  I=N+]-K  which  increases 
from  1-1  at  the  towpoint  to  l-N  at  the  towed  body.  The  first  three  data  cards 

differ  from  those  of  OMWAY  and  specify: 

1  ' 

1.  Name  of  program,  OMFLO. 

2.  Steady-state  body  lift  and  drag  BLO  and  BDO,  and  body  weight  BW. 

3.  Magnitude  and  phase  of  oscillating  lift  and  drag  forces  on  the  body,  BLM , 
DBL,  DBM,  DBD,  and  frequency  OM. 

These  data  are  written  out.  Do  is  computed  and  Equation  [14]  is  used  to  compute 

X,  which  is  represented  by  AMC/ULM.  Equations  [33]  and  [34]  are  used  to  find 

' 
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T0  and  <Po  at  K=1  and  Vo,  Vo,  0o.  To,  xo,  and  yo  are  computed  by  using  Equations 
[12],  [13],  [19],  [20],  [31],  [32],  [SS],  and  [57]  and  are  written  out  as  functions 
of  J.  0 o  and  To  are  relisted  as  functions  of  7,0/^  is  set  real  and  Tja,  Via,  *nd  Vja 
are  set  equal  to  zero  at  7s/,  and  the  A-mode  solution  is  computed  by  integrating 
Equations  [35]—  [38]  on  1  <I<>N.  Then  Tjb  is  set  real  and  0/a,  V\b,  and  Vjb  are 
are  set  equal  to  zero  at  7s  1  and  the  5-mode  solution  is  computed  in  the  same  way. 
Next  00,  To,  Via,  Vja>  4ia>  Tia,  Vjb,  V\b,  0/*,  and  Tib  are  re-listed  as 
functions  of  K  instead  of  7.  The  required  values  of  P  and  Q  at  Af=l  are  found  from 
Equations  [39]  and  [40]  and  Ca  and  Cb  are  found  from  Equations  [51]  and  [52]. 

V ),  V i ,  0/ ,  and  Tj  are  found  along  the  cable  from  Equations  [45] -[48]  and  their 
magnitudes  and  phases  are  computed  and  written  out.  Finally  Equations  [56]  and 
[58]  are  used  to  And  xj  and  yj,  and  their  magnitudes  and  phases  are  computed  and 
written  out. 

Representative  times  for  a  run  involving  a  cable  nine  integration  steps  long  (of 
ten  feet  each)  are  28  seconds  compilation  time  and  14  seconds  computation  time  on 
a  CDC  6700  computer. 

Certain  modiAcations  to  the  programs  are  easily  made.  For  example,  the 
programs  listed  require  angles  in  radians  and  any  consistent  mass-length-time  units  for 
dimensional  quantities.  They  have  been  adapted  to  treat,  as  input  and  printout 
quantities,  angles  in  degrees,  masses  in  pounds,  lengths  in  feet,  and  speeds  in  knots. 
They  have  also  been  adapted  to  examine  behavior  over  a  range  of  frequencies,  by 
repeating  the  unsteady  part  of  the  program  in  a  DO  loop,  with  lowest,  highest,  and 
incremental  frequencies  speciAed  on  an  additional  data  card. 

A  computer  experiment  was  performed  on  the  Quadrant  1  program,  OMWAY, 
by  examining  the  frequency  range  0.01  to  0.80  hertz  in  steps  of  0.01  hertz  for 
cable  lengths  of  200,  400,  600,  800,  1000,  and  1200  feet  at  tow  speeds  of  6,  10,  and 
14  knots.  As  criteria  for  successful  performance  it  was  required  that  the  perturbation 
quantities  must  remain  below  the  steady-state  quantities.  It  was  found  that  Vm,  the 
magnitude  of  the  tangential  velocity,  at  (or  very  near)  the  towpoint  was  the  quantity 
which  failed  at  the  lowest  frequency.  The  only  cable  length  where  its  behavior  was 
satisfactory  throughout  the  frequency  and  speed  ranges  was  200  feet;  all  the  longer 
ones  failed  but  showed  better  behavior  with  increasing  speed.  The  600-foot  cable 
failed  at  0.26,  0.34,  and  0.48  hertz  at  6,  10,  and  14  knots  while  the  1200-foot 
cable  failed  at  0.13,  0.20,  and  0.30  hertz  at  the  same  speeds.  Values  for  the  8 00- and 
1000-foot  cables  lie  in  between  these  and  exhibit  the  same  monotonic  relation  to  tow 
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m**1 


, _ • 


aa 


speed.  The  400-foot  cable  failed  at  0.49  and  0.57  hertz  at  6  and  10  knots  and 
behaved  satisfactorily  at  14  knots. 

The  next  most  sensitive  quantity  was  the  perturbation  tension  at  or  very  near 
the  ship.  It  behaved  well  throughout  the  frequency  and  speed  ranges  for  cable 
lengths  of  800  feet  and  less.  For  the  1200  foot  cable  it  failed  at  0.46,  0.46,  and 
0.61  hertz  at  6,  10,  and  14  knots. 

The  qualitative  conclusions  are  that  the  small-perturbation  frequency-domain 
analysis  is  valid  for  short  cables  undergoing  low-frequency  oscilations  and  that  increas¬ 
ing  towspeed  can  have  a  stabilizing  effect.  Apparently  as  the  forcing  frequency  ap¬ 
proached  0.80  hertz,  a  resonant  region  was  being  approached;  frequencies  high 
enough  to  be  past  such  a  region,  because  the  inertia  of  the  cable  and  body  prevent 
excitation  of  their  motion,  might  be  beyond  the  range  of  practical  interest. 
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TOWED  BODY 

Figure  1  -  Quadrant  I  Towing  Configuration 


TOWING  SHIP 


Figure  2  -  Quadrant  II  Towing  Configuration 
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APPENDIX  A 

DERIVATION  OF  EQUATIONS  OF  MOTION 

Equations  [1]— [4],  the  equations  of  cable  motion,  can  be  derived  by  considering  a 
small  element  of  the  cable  shown  in  Figure  1  or  Figure  2.  If  x  and  y  represent  the 
coordinates  of  a  point  on  the  cable  then  x  and  y  are  functions  of  s  and  t  which 
satisfy 


3*  . 

—-*Cos0 

bs 


I A 1  ] 


by 

■57-“"*  (A2I 


bx  by 

—  sin0-““  cos0  ■  U 
bt  bt 


[A3] 


and 


—*  cos0  +  -7-  sin0  *  V 
bt  bt 


[A4] 


Take  the  derivative  of  Equation  [A3]  with  respect  to  s,  interchange  the  order  of  the 
s -  and  t-differentiations  of  x  and  y,  and  substitute  from  Equations  [Al]  and  [A2]  to 
obtain 

bU  3  ,  bx  30  3  ...  by  30 

T~  =  ~  (cos 0)  sin 0  +  —  cos0  ~~  ~~  (sm0)  cos0  +  —  sin0  — 

bs  bt  bt  bs  bt  bt  bs 


which  becomes  Equation  [1]  after  substitution  of  Equation  [A4].  In  the  same  way, 
take  the  derivative  of  Equation  [A4]  with  respect  to  s,  interchange  the  order  of 
differentiations,  and  substitute  to  obtain 

bV  b  ,  bx  .  b<f>  b  by  30 

~  =  77  (cos0)  cos0  -  —  sm0  —  +  —  (sin0)  sin0  +  -7  cos0  — - 

bs  bt  bt  bs  bt  bt  bs 

This,  with  the  substitution  of  Equation  [A3],  is  Equation  [2]. 

The  apparent-momentum  vector  of  a  segment  of  length  6s  has  components 
p6s(l+X)C/,  n6 sV,  0  in  a  right-handed  orthogonal  coordinate  system  with  axes  in  the 
directions  of  U,  V,  and  upwards  from  the  paper  of  Figure  1  or  Figure  2.  Here 
M  is  the  mass  of  the  cable  per  unit  length  and  pX  is  the  added  mass  in  water  for 
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motion  in  the  direction  of  U.  The  angular  velocity  has  components  0,0,  bffbt  in  the 

. 

fixed  coordinate  system  defined  by  the  x  and  y  axes  and  the  direction  up  out  of  the 
paper  in  either  figure.  Therefore  the  rate  of  change  of  apparent  momentum  observed  from  the 
fixed  coordinate  system  in  either  figure  has  components  along  the  directions  of  V  and  V  equal 
to 


bU  30 

M  Ml  +X)  —  — 


and 


bV  30 

+n6s(\  +X){/  — 


Therefore  the  dynamic  equations  of  motion  of  the  cable  segment  are 


[3  U  30  1 

(,+»F-KrJ 

[bV  30 1 

57  ♦<l+x,£'irJ 


=  resultant  force  on  6s  in  the  direction  of  U  [AS] 


=  resultant  force  on  6s  in  the  direction  of  V  [A6] 


The  right-hand  sides  of  Equations  [AS]  and  [A6]  are  easily  computed  and  the  results 
are,  after  division  by  6s,  Equations  [3]  and  [4]. 
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APPENDIX  B 

LISTING  AND  SAMPLE  OUTPUT  OF  OMWAY 


5 


1° 


IS 


20 


25 


30 


35 


40 


45 


50 


55 


PROGRAM  OMWAY ( INPUT .OUTPUT . TAPES* INPUT  »TAPE6*OUTPUT ) 

DIMENSION  PHO(400) »  T0(400).  DO (400).  VO (400).  PHRAUOO). 

♦  P HRS (400).  PHIAI400) .  PHIB(400)»  TRA(400)»  TRBI400).  TIA(400)» 

♦  TIB (400) t  URA(400)t  URB(400).  UlA(400).  UIBI400).  VRA(400). 

♦  VRB(400).  VIA(400) .  VIB(400) .PHR(400) *  PHK400).  TR(400). 

♦  TI(400).  UR (400).  U1 (400) *  VH(400).  VI(400).  FL<400).  F6(400). 

♦  DFL (400) »  OFG(400) .  T(B) 

COMMON  DS.WUL.D. ALAMO. ALAMl  .ALAM2.BLAM1 . 

1BLAM2. AGAMO . A6AM1 .AGAM2.BG AMI .BGAM2. 

2OM.ULM.APP.C.PH0.OFL.DFG. J.TO.PHR.PHI .FL.FG 
EXTERNAL  OAUX.OARN 
RE AO  (5.11)  TITLE. 

♦  BL  *  BO  *  BW. 

♦  ASM.  OA$.  AHM.  OAH.  OM. 

♦  BM.  BMS.  BMH. 

♦  WUL*  ULM.  OS.  N. 

X  ♦  C.  ON.  CM.  CO.  AMP. 

♦  ALAMO.  ALAMl.  ALAM2.  BLAM1.  BLAM2. 

♦  AGAMO.  AGAM1.  AGAM2.  BGAMI.  BGAM2 

C  WRITE  TITLE  ANO  INPUT  DATA 

WRITE  (6.21)  TITLE 

WRITE  (6.22)  BL  .  BO  .  BW 

WRITE  (6.23)  ASM.  DAS.  AHM.  DAH.  OM 

WRITE  (6.24)  BM.  BMS.  BMH 

WRITE  (6.25)  WUL.  ULM.  ,DS.  N 

WRITE  (6.26)  C.  ON.  CH.  CO.  AMP 

WRITE  (6.27)  ALAMO.  ALAMl.  ALAM2.  BLAM1.  BLAM2 

WRITE  16.28)  AGAMO.  AGAMl.  AGAM2.  BGAMI*  BGAM2 

C  INITIALIZE  , 

D*  0»5*DN*C*C«CH#CO 

AMC •  3. 1415R3*CH*CM*DN* AMP/4.0 

APP  ■  1.0  ♦  AMC/ULM 

CXO  ■  0.0 

CYO  *  0.0 

CXR  ■  0.0 

CXI  *  0.0 

CYR  ■  0.0 

CYI  *  0.0 

JA  ■  N-l 

C  SOLVE  STEADY-STATE  PROBLEM 
ERR*. 00001 
ERA*. 001 
CALL  ATA(BO.YK) 

PHOB  ■  AT AN2 ( BM-BL . YK  > 

TOB  *  SORT (80*80* (BW-BL) •(BW-BL> ) 

WRITE  (6.31) 

WRITE  (6.32) 

PH0(1)  *  PHOB 
10(1)  *  TOB 

C  MULTIPLY  TO  FINO  VELOCITY  COMPONENTS 
DO  101  J  ■  l.N 
U0(J)  *  C*5IN(PH0 (J) ) 

V0(J>  ■  C*COS(PHO(J) ) 

WRITE  (6.33)  PHO(J).  TO(J).  UO(J).  V3(J).  CXO.  CYO.  J 
IF  (J.EQ.N)  GO  TO  101 


25 


60 


6S 


70 


75 


SO 


85 


90 


95 


100 


105 


110 


C  INTEGRATE  TO  FIND  TENSION*  ANGLE*  ANO  CABLE  CONFIGURATION 
T(1)*T0(J> 

T(2)«PH0(J) 

TO)  >  CXO 
T (A)  «  CYO 

FIRST*. 0 
ALS*FLOAT  (J)*DS 

CALL  KUTMEH(4.ALS*T»ERR»DS»FIRSI»HCX»ERA»DAUX> 

T0(J*1)«T(1> 

PH0(JM)»T(2> 

CXO  •  T  (3) 

CYO  ■  T (4) 

101  CONTINUE 
C  SOLVE  DYNAMIC  PROBLEM 
WRITE  (6*34) 

WRITE  (6*35) 

C  COMPUTE  VISCOUS-FORCE  TERMS  ALONG  CABLE 
00  201  J  ■  1*  JA 

FLU)  *  D*  ( ♦  AL  AMORAL  AMI  *COS  (PMO  U) )  ♦ALAM2*C0S (2»0*PH0 ( J) )  ♦ 

♦  BLAM1*SIN(PH0 ( J) ) ♦BLAM2*SIN(2.U*PH0 (J) ) ) 

F6( J)  -  0*1  AGAMO*AGAMl*COS (PMO ( J) ) ♦AGAM2*COS(2.0*PH0 ( J) ) 

♦  *BGAM1*SIN(PH0 (J) ) ♦BGAM2*SIN(2.0*PH0 ( J) ) ) 

OFL(J)  *  D*(-ALAM1*SIN(PH0(J) )-i£»0*ALAM2*SIN(2.0*PHO(J) ) 

♦  ♦BLAM1*C0S (PHO (J) ) ♦2.0*BLAM2*CUS(2.0*PH0 ( J) ) ) 

DFG(J)  *  U* (-AGAM1*SIN (PHO ( J) ) -2»0*AGAM2#SIN(2.0*PH0 ( J) ) 

♦  *  BGAM 1*COS(PHO(J) )*2.0 *BGAM2*CUS ( 2 . 0  *PHO ( J  > ) > 

201  CONTINUE 

C  COMPUTE  MODAL  SOLUTIONS 
C  COMPUTE  A  -MODE  SOLUTION 

C  ASSIGN  VALUES  TO  VELOCITY  COMPONENTS  AT  BOTTOM 
URA(l)  *  0.000001 
U1A(1)  *  0.0 
VRA(l)  a  0.0 
VIA(l)  «  0.0 

C  COMPUTE  CORRESPONDING  VALUES  OF  ANGLE  ANO  TENSION  AT  BOTTOM 
PHRA(l)  a  OM* (-UI A ( 1 ) * (8M+BMS*SIN(PH0 ( 1 ) )*SIN(PHO (1) ) ♦ 

♦  BMH*COS (PHO (1>) *COS (PHO  ( 1 > ) >  «V1 A ( 1 ) • (BMH-BMS) *SIN (PHO ( 1 > ) *COS ( 

♦  PHO (1)))/I0(1> 

PHIA(l)  a  -OM*(-URA(1)*(BM*BMS*SIN(PHO(1) ) *SIN(PHO (1 )) ♦ 

♦  BMH*COS (PHO ( 1 ) ) •COS(PHO  ( 1 ) ) ) ♦ VKA ( 1 ) • (BMH-BMS) *SIN  (PHO ( 1 ) ) «COS ( 

♦  PHO ( 1 > ) ) /TO ( 1 ) 

TRA(l)  «  OM*(-UIA(1)*(BMH-BMS)«SIN(PHO(1))*COS(PHO(1)) 

♦  *VI A ( 1 ) * (BM*BMS*COS (PHO ( 1 ) ) *COS (PHO ( 1 ) ) ♦ 

♦  BMH*SIN (PHO ( 1 ) ) *SIN(PHO (1 ) ) ) ) 

TIA(l)  a  -OH* (-URA ( 1 )* (BMH-BMS) *SIN(PHO ( 1 ) ) *COS(PHO (1 ) ) 

♦  *VRA ( 1 ) * (BM*BMS*COS (PHO  ( 1 ) ) *COS (PHO ( 1 ) ) ♦ 

♦  BMH*SIN (PHO ( 1 ) )*S1N(PH0(1) ) ) ) 

C  FIND  VELOCITY  COMPONENTS*  ANGLE*  AND  TENSION  BY  INTEGRATING  UP 
C  THE  CABLE 

DO  202  J  «  1*  JA 
T(1)«PHRA(J) 

T(2)aPHIA(J) 

TO)aTRA(J) 

T(4)»TIA(J) 

T(5)>URA(J) 
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115 


120 


125 


130 


135 


140 


145 


150 


155 


160 


165 


T (6)*UIA( J) 

T(7)*VRA(J> 

mxviMji 

FIRST*. 0 
ALS*FLOAT ( J)*DS 

CALL  KUTMER(8*ALS*T »ERR«0S»FIRS1  *HCX* ERA* DARN) 

phra(j*d*t<i) 

PHIA( J*1)b1 (21 
TRA<J*l)*T(3) 

TIA(J*1)*T (41 
URA(J*1)*T(5) 

UIA(Jd)*T<6> 

VRA(J*1)«T(7) 

VI  A<  J«-l  >  «T  «8I 
202  CONTINUE 

C  COMPUTE  B  -MODE  SOLUTION 

C  ASSIGN  VALUES  TO  VELOCITY  COMPONENTS  AT  BOTTOM 
'  URO(l)  -  0.0 
UIB(l)  ■  0.0 
VRB(l)  ■  0.000001 
VIB(1)  *  0.0 

C  COMPUTE  CORRESPONDING  VALUES  OF  ANGLE  ANO  TENSION  AT  BOTTOM 
PHRB ( 1 )  ■  OM*(-UIB(l>*(BM*BMS*SIN(PHO(n)*SIN(PHO(l))« 

♦  BMH*COS <PM0  ( 1 ) » <*COS  (PHO  ( 1 ) ) )  ♦V IB  ( 1 )  •  (BMH-BMS )  *SIN IPMO  ( 1 )  > «*COS  ( 

♦  PH0(1) ) )/T0(l) 

PHIB(1 )  ■  -OM*(-URB(  1>  #(BM*BMS*SIN(PH0(1)  )*SIN(PHO  (1 ) )  ♦ 

♦  BMH*COS  (PMO  (1)1  *COS (PHO  (111)  ♦  VKB  ( 1 » •  (BMH-BMS)  *SIN  (PHO  (111  *COS  ( 

♦  PHO (11)1/10(1) 

TRB(l)  -  OM*(«UIB(l)*(BMH— BMS)#SIN(PM0(1)  )*COS(PHO(1) ) 

♦  ♦  VIB  ( 1 1  •  (HM*BMS*COS  (PHO  ( 1 )  )  *C05»  (PHO  ( 1 ) )  ♦ 

♦  BMH*SIN(PHO(1))*SIN(PHO(1)))) 

TIB(l)  ■  -OM*(-UR0(l)*(BMH-BMS)*SIN(PHO(l))»COS(PHO(m 

♦  ♦  VRB  ( 1 )  *  ( BM*HMS*COS  ( PHO  ( 1 )  ) *COS (PHO  ( 1 ) )  ♦ 

♦  BMH*S1N(PH0(1) )  *SIN  (PHO (1) ) ) ) 

C  FIND  VELOCITY  COMPONENTS*  ANGLE*  ANO  TENSION  BY  INTEGRATING  UP 
C  THE  CABLE 

DO  203  J  ■  1*  JA 
T(1)*PHR8(J) 

T (2)*PHIB( J) 

T(3)*TRB(J) 

T(4)»TIB(JI 

T(5)»URB(J) 

T(6)*UIB(J> 

T(7)«VRH(J) 

T (8)*V18( J) 

FIRST*. 0 
ALS*FLOATU)«OS 

CALL  KUTMEH(B«ALS*T*ERR*OS*FlKST*HCX*ERA*OARN> 

PMRB(J*1)*T(1) 

PMIB(J*1)*1 (2) 

TRB(J*1)*T (3) 

TIB(J»1)*T(4> 

UR8(JM>*T(5> 

UIB( J*1 >*T (6) 

VRB  ( *)♦!  )*T  (7) 
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VIB ( J*l ) «1 (0) 

203  CONTINUE 

C  SUPERIMPOSE  MODAL  SOLUTIONS 
C  COMPUTE  VELOCITY  COMPONENTS  FORCEO  AT  TOP 
CSR*OM«ASM«SIN(OAS) 

CSI»-OM«ASM«COS(DAS) 

CHR«OM#AHM*SIN(DAH) 

CHI «-OM*AMM«COS(OAH) 

UR (N)  ■  CS«*SIN(PHO(N) T  -  CMR*COS(PHO (N) T 

Ul  (N)  =  CSI*SIN(PHO (N) )  -  CHI*COS(PHO(N>> 

VR<N)  =  CSR*COS (PHO (NT )  ♦  CHR*SIN(PHO (N) ) 

VI  (N)  «  CSI*COS(PHO(N) )  ♦  CM I •SIN (PHO (N) ) 

C  COMPUTE  CONSTANTS  NEEDED  TO  SUPERIMPOSE  MODAL  SOLUTIONS 

DENOMK  *  URA (NT *VRB (N)  -  UIA (N) *VIB(N)  -  URB(N)*VRA(N)  ♦ 

♦  UIB (N) *V I A (N) 

OENOMT  >  UIA (NT *VRB  < N >  ♦  URA (N) •VIB(N)  -  UIB (N) *VRA (N)  - 

♦  URB (N) *VIA (N) 

ANUMK  *  UR  (NT  *VRB(N)  -  UI  (N)*VI0(N)  -  URB(N)*VR  (N)  ♦ 

♦  UIB(N)*VI  (NT 

ANUMI  s  ui  (N)*VRB(N)  ♦  UR  (N)*VIB(NT  -  UIB(N>*VR  (NT  - 

♦  URB (NT  *VI  (N) 

BNUMR  *  URA (N) *VR  (N)  -  U1A(N)*VI  (NT  -  UR  (N»*VRA(N)  ♦ 

♦  UI  (N)*V1A(N) 

BNUMI  »  UIA(N)*VR  (NT  ♦  URA(N>«VI  (NT  -  UI  (N)*VRA(NT  - 

♦  UR  (NT  *VIA(NT 

DENOM2  ■  DENOMR*DENOMR  ♦  DENOMl*DENOMI 
CAR  *  ( ANUMR*l)ENOMR  ♦  ANUMI  *DENOMI )  /0EN0M2 

CAI  >  (-ANUMR*OENOMI  ♦  AN(/MI*DEN0MR)/DEN0M2 
CBR  «  (BNUMR*DENOMR  ♦  BNUMI *DEN0MIT/DEN0M2 
CBI  *  ( -BNUMR*DENOM I  ♦  BNUM I *DEN0MR I /DEN0M2 
C  COMPUTE  VELOCITY  COMPONENTS*  ANGLE*  AND' TENSION  ON  CABLE 
DO  204  J  «  1  *N 

UR ( J>  «  CAR*URA ( JT  -  CAI*UIA(J)  ♦  CBR*URB(JT  ■-  CBI*UIB(JT 

UI(JT  «  CAI*URA ( J)  ♦  CAR*UI A ( J)  ♦  CBI»URB(JT  ♦  CBR*UIB(J» 

VR( JT  ■  CAR*VHA ( J)  -  CAI*VIA(JT  ♦  CBR*VRB(JT  -  CBI»VIB(JT 

VI (JT  «  CAI*VHA ( JT  *  CAR*VI A ( J)  ♦  C0I»VRB(JT  ♦  CBI*VIB(JT 

PHR(J)  *  CAR*PHRA ( JT  -  CAI*PH1A(J)  «  CBR*PHRB(J>  -  CBI*PHIB(J) 
PHI (J)  ■  CAI*PHRA ( J)  *  CAR*PH1A( JT  *  CBI*PHRB(JT  ♦  CBR*PHIB(J) 
TR( JT  «  CAR*TRA(J»  -  CAI*TIA(J)  ♦  CBR*TRB( JT  -  CBI*TIB<J) 

TI(JT  >  CAI*TRA ( J)  ♦  CAR*TIA(J>  *  CBI*TRB(JT  ♦  CBR*TIB(JT 

C  COMPUTE  PHYSICAL  VELOCITY  COMPONENTS  ON  CABLE 
URP  «  C*PHR(JT*COS(PHO(JT  T  ♦  UR(J! 

UIP  =  C*PHI  ( J)  *COS  (PHO  ( JT  >  ♦  UKJJ 
VRP  «  -  C*PHR(JT*SIN(HHO(JT  T  ♦  VR(J) 

VIP  *  -  C»PHI(JT*SIN(PHO(JTT  ♦  VI(JT 

C  L?*?UTE  MAGNITUDE  AND  PHASE  OF  VELOCITY  COMPONENTS*  ANGLE*  AND 
C  TENSION  ON  CABLE 

UM  ■  SORT (URP*URP  ♦  UIP*U1PT 
CALL  ATA(URP*  YK) 

OU  «  ATAN2 (UIP*  YKT 

VM  »  SORT  (VRP* VRP  ♦  VIP*VIP! 

CALL  AT A (VRP*  YK) 

OV  *  ATAN2 ( VIP*  YK) 

PHM  «  SORT (PHR(JT*PHR( JT  ♦  PHI ( JT *PHI ( JT ) 

CALL  ATA(PHRIJ)*  YK) 
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225 

230 

235 

2*0 

245 

250 

255 

260 

265 

270 

275 


DPH  >  ATAN24PHI ( J) *  ii w 

TM  ■  SORT (TR4JI*TR4J)  ♦  TI4J)*TIIJ)) 

CALL  ATA4TR4J) ♦  YK) 

DT  ■  ATAN24TI4J) *  YK) 

WRITE  (6*36)  PHM.  DPH.  TM*  DT*  UN*  DU.  VM»  DV*  J 

204  CONTINUE 
WRITE  46*37) 

C  COMPUTE  MAGNITUDE  A  NO  PHASE  OF  DYNAMIC  CABLE  CONFIGURATION 
00  205  J  ■  1*  N 
CXM  >  SORT (CXR*CXR*CXI*CXI ) 

CALL  ATA4CXR*  YK) 

OCX  •  ATAN24CXI*  YK) 

CYN  a  SORT <CYR*CYR*CYI*CYI) 

CALL  ATA4CYR*  YK) 

OCY  a  ATAN24CYI *  YK) 

WRITE  46*36)  CXM*  OCX*  CYM*  OCY*  J 
IF  4J.EQ.  N)  60  TO  205 
CXR  a  .  PHR4J)»SIN4PHO(J))*OS 
CXI  a  .  PHI ( J)*SIN(PH0( J) )*OS 
CYR  a  PHR(J)*COS4PH04 J) )*DS 
CYI  >  PHI  4 J) *COS 4PH0  4 J) ) #OS 

205  CONTINUE 

II  FORMAT  4  10X*  A6  /  I0X*  3F10.5/  10X*  5F10.5/  10X*  3F10.S/ 

♦  10X*  3F10.5*  II0/I0X.5F10.5  /  10X*  5FI0.5/  10X*  5F10.5) 

21  FORMAT  41H1*  10X*  57HSURFACE  SHIP  IN  A  SEAWAY  TOWING  A  DEEPLY 

♦  SUBMERGED  WEIGHT/  1H-.  lOX*  A6///I 

22  FORMAT  4  IX*  8HBL  ■  «  FI0.5.  4X*  8MBD  a  .  F10.5*  4X* 

♦  SMBW  a  ,  F10. 5  //) 

23  FORMAT  4  IX*  AHASM  a  ,  F10.5*  4X»  8H0AS  a  ,  FI0.5*  4X* 

♦  BHAHM  a  ,  F10.5*  4X«  8HOAH  a  ,  F10.5.  4X*  8M0M  a  , 

♦  F10. 5  //) 

24  FORMAT  4  IX*  8HBM  a  ,  FiO.S*  4X*  8HBMS  a  ,  F10.5*  4X« 

♦  BHBMH  a  ,  FIO.S  //) 

25  FORMAT  4  IX.  BHMUL  a  ,  F10.5*  4X*  8HULM  a  *  FIO.S*  4X* 

♦  8HDS  a  .  FIO.S*  4X»  8MN  ■  .  110  //) 

26  FORMAT  4  IX.  SMC  a  ,  FIO.S*  4X*  8H0N  a  *  FIO.S.  4X* 

♦  AHCH  a  .FIO.S.  4X*  8HCO  a  .FIO.S.  4X.  8HAMP  a  .FIO.S//) 

27  FORMAT  (  IX.  8HALAM0  a  .FIO.S*  4X*  BHALAM1  a  .  F10.5*4X« 

«  8HALAM2  a  ,  FIO.S*  4X*  8HBLAM1  a  .  FIO.S*  4X*  8HBLAM2  a  , 

♦  FIO.S  //) 

28  FORMAT  (  IX.  8HAGAM0  a  .  F10.5.  4X*  BHAGAM1  a  *  FIO.S*  4X* 

♦  BHAGAM2  a  ,  FIO.S*  4X*  8HBGAM1  a  ,  F10.5*  4X»  8HBGAM2  a  . 

♦  F10.5  //) 

31  FORMAT  '  22H1  STEADY-STATE  SOLUTION) 

32  FORMAT  (  1H0*  IX.  6HPH0 ( J)  *  10X*  6HT04J)  .  10X*  6HU04J)  *  10X* 

♦  6HV04J)  .  10X.  6HCA0  ( J)  *  10X*  6MCY0  ( J)  *  10X*  1HJ  //) 

33  FORMAT  4  2X.  F9.6*  7X.  F13.6*  3X*  24  F11.6*  57  .  24F13.6*  3X>  *  14) 

34  FORMAT  4  17H10YNAM1C  SOLUTION) 

35  FORMAT  4  1H0.  IX*  6HPHMU)*  9X»  6HDPH4J).  9X.  6HTM4J)  *  9X* 

♦  SHOT 4 J)  .  9X.  6HUM4J)  .  9X*  6HUU4J)  .  9X*  6HVM4J)  .  9X*  6HOV4J)  * 

♦  9X.  1HJ  //) 

36  FORMAT  4  2X.  24  F9.6*  6X).  F12.6*  3X*  F9.6*  6X.  24F10.6*  5X* 

♦  F9.6*  6X).  14) 

37  FORMAT  4  1H1.  IX*  6HCXM4J) .  10X*  6HOCX4  J) .  10X*  6MCYM4J).  I0X* 

♦  6HDCY4J).  10X.  1HJ  //) 

38  FORMAT  4  2X.  24  F12.6*  4X*  F9.6*  7X>*  14) 

STOP 

END 
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SUBROUTINE  ATA(X*YK) 

C  THIS  SUBROUTINE  PREVENTS  ATAN2  f ROM  DIVIDING  BY  ZERO  WHEN 
C  COMPUTING  THE  ARCTANGENTS  Of  PI/2  AND  3*RI/2 
ER  >  0.00001 

If  jtABS(X).GT.ER)  GO  TO  10 

YK  a  £H  l 

RETURN 

10  YK  a  X  '  •’  ' 

RETURN  i  ’  ' 
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20 


25 


30 


i 


35 


AO 


AS 


50 


55  ' 


/ 


SUBROUTINE  KUTMER(N*T*V0*EPS*H*f IRST.HCX.A.DAUX) 

C  THIS, SUBROUTINE  PERfORMS  fOURTH-OHOER  XUTTA-MERSON  INTEGRATIONS 
.  DIMENSION  YO ( 10) « Y1 110) *Y2< 10) tf 0 (10) *fl (10) *P2( 10) 

If (flRST|20. 10.20 
10  HC-H 
IPLOC-1 

fIRSTal.  I  1  I  I 

20  LOCaO 

.  HCX*HC  '  i  l 

30  CALL  OAUX(l*YO*fO)  '  , 

39  DO  AO  lal.N  ' 

AO  Ylil)aVO(I)«(HC/3.)*fO(I) 

CALL  DAUX (T*HC/3. * VI *fl )  ,  i 

DO  SO  lal.N 

50  Ylll)aYO(I) ♦ (HC/6. ) *f 0 ( I ) ♦ (HC/6. ) *fl 1 1 ) 

CALL  DAUX (T*HC/3. » Y1 *fl )  / 

00  60  lal.N  ' 

60  Y1(I)«Y0(I) ♦HC/B.#f 0 (I ) ♦•375#HC*f 1 (I )  / 

,  CALL  DAUX n*HC/2. • VI »f2l 
00  70  lal.N, 

70  Vl(I>aY0(I)*HC/2.*f0(l)-1.5*HC*H(I>*2.*HC*f2(I> 

CALL  OAUX11  «HC».Y1  *fl )  1 

'  DO  60  lal.N  , 

80  Y2( I )aY0 (I) ♦HC/6.*f 0 1 1) ♦•666666b7*HC*f 211) ♦ (HC/6. )*f 1 (I) 

INC-0 

DO  110  I°1.N  i  I 

ZZZaABSlYl (I) )-A 

1  If(ZZZ)  85.A7.S7  1 

85  ERROR* ABS ( .2*  ( Y 1  ( I )  -Y2  ( I ) ) ) 

If (ERROR-A) 100.100*90  i 

87  ERRORaABS ( • 2“*  2* Y2 1 1 ) /Y 111)) 

.  If (ERROR-EPS) 100.100*90  .  . 

90  X«128.aABS(HC)-ABS(H)  1 

If(X)  91*95*95 

91  WRITE (6*92)  T. ERROR 

92  f ORMAT (21H  RELATIVE  ERROR  AT  Ta  1P1E12.3.3H1S  fl0.6  /  > 

'  fIRST  a  2.  1  i  l 

RETURN 
95 . HC“HC/2. 

IPL0Ca2  *IPL0C  ' 

L0Ca2  *L0C 
HCXaHC 
GO  TO  30 

100  If (ERR0R*6A.-EPS) 110*110*101 

101  INCal  ,  i  1 

110  CONTINUE  1 

111  T«T*HC 

DO  112  I>1*N  i 

112  Y0(I>-V2(I> 

LOC»LOC*l 

If 1LOC-IPLOC) 120*210*210  1  , 

120  If (INC) 210 *130 *210  ' 

’130  If (L0C-1L0C/2) *2)210. 1A0. 210 
1A0  If (IPLOC-D210.210.200  / 

200  HCa2.*HC 
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1  KUTM0050 
KUTH0060 
KUTM0070 
KUTM0080 
KUTH0090 
KUTM0100 
KUTM01 10 
KUTN0120 
KUTM0130 
KUTMOIAO 
KUTM0150 
KUTM0160 
KUTM0170 
KUTM0180. 
KUTH0190 
KUTM0200 
KUTM0210 
KUTH0220 
KUTH023Q 
KUTH02A0 
KUTM0250 
KUTM0260 
KUTM0270 
KUTN0280 
KUTM0290 
KUTH0300 
KUTM0310 
KUTH0320 
KUTM0330 
KUTM03A0 
KUTM0350, 
KUTM0360 
KUTM0370 
KUTM0380 
KUTH0390 
KUTM0A00 
KUTH0A10 
KUTM0A20 
KUTM0A30 
KUTMOAAO 
KUTM0A50 
KUTM0A60 
KUTM0A70 
KUTM0A80 
1  KUTM0A90 
KUTM0500 
KUTM0510 
KUTH0520 
KUTM0530 
KUTM05A0 
KUTH0550 
( 


I 


I 


I 


l 


LOCaLOC/2 

IPL0Ca!Pl0C/2 

210  IF (IPL0C-L0C >30*220*30 
220  RETURN 

60  i  END 


SUBROUTINE  DAUX(Z*TX*F) 

C  THIS  SUBROUTINE  PROVIDES  THE  INTEGRANOS  FOR  KUTMER  WHEN  IT 
C  COMPUTES  THE  STEADY-STATE  TENSION*  ANGLE*  AND  CABLE  CONFIGURATION 
i  DIMENSION  TX (8) *  F(8I *PH0(400) *DFL(400) .OFGUOO) *T0 (600) * 

5  ♦  PHRI200) «  PHI (200) 

COMMON  OS*MUL*D*ALAMO*ALAMI *ALAM2*BLAH1* 

I BLAM2 * AGAMO  *  AG AMI • AGAM2 • BGAM I  * BGAM2 * 

I  20M*ULM*APP*C*PH0*DFL*DFG«J*T0*PhR«PHI 

TXl-TXU) 

TX2«TX(2) 

C1>C0S(TX2) 

S1»SIN(TX2) 

C2«C0S<2.*TX2) 

S2>SIN(2.*1X2) 

FL*D# (ALAMO *AL AMI *C1* ALAM2#C2*BL AMI *S1 *BLAM2*S2) 

FG>  D* ( AGAMO* AG AM 1  *C  1  ♦ AGAM2*C2*BGAMI*S I *BGAM2*S2 1 
IFITX1.EQ.0.)  GO  TO  SO 
F (2)*  (WUL*C1-FL)/TX1 
F(3)  ■  Cl  ' 

F(4)  «  ,  SI 

F  (1 ) ■  (WUL*S1*FG) 

RETURN 
F (2)*0.0 
GO  TO  7 
END 


IS 

I 

! 

20 

7 

50 

i 

25 


KUTM0560 

KUTM0S70 

KUTM0540 

KUTM0S90 

KUTM0600 


i 


1? 


IS 


20 


,  i 

25 


30 


SUBROUTINE.  DARN(Z.TtF) 

C  THIS  SUBROUTINE  PROVIDES  THE  INTEGRANOS  FOR  KUTMER  WHEN  IT 
C  COMPUTES  THE  DYNAMIC  ANGLE*  TENSION*  AND  VELOCITY  COMPONENTS 
DIMENSION  T  (B)  *  FIB)  »PH0  (400)  tUFLUOO)  *DFG(400)  *TO  (400)  * 

♦  FL(AOO) *F6(400) *PHR(600) »PHI (400) 

COMMON  OS*MUL*0« ALAMO *ALAM1 « ALAM2«BLAM1 • 

1 BL  AM2  «  AG AMO • AO AM 1  *  AG AM2  «  BG AM 1 «  BGAK2  * 

20M*ULM«APP*C*PHO»DFL»DFG»  J*TO»PHR»PHl  »FL*FG 
C)*COS(PHO( J) ) 

S1«SIN(PH0 ( J) ) 

IF  (TO  (J)  .Ell. 0.0)  GO  TO  50 

F(1)«M-OM*ULM*(APP*T(6)*(APP-1.V)»C»T!2)pCI)-T(3)<MPHOU*I)-PMO(J) 

♦)/DS-(T(l)*(T(5)*Cl-T(7)*Sl)/C)*DFL(J) 

♦  -2*0*FL( J)*(T (5)*S1*T (7)*CI)/C-WUL*T (1)*SI)/T0IJ) 
F<2)«(  0M*ULM* (APP*T (5) ♦ ( APP-I «0) *C*T ( I ) *Cl ) -T (*)*<PHO ( J*l) -PHO ( J) 

♦)/DS-(T (2) ♦ (T (6)*C1-T(8) *S1 )/C)#DFL(J) 

♦  -2.0'tFL(J)»(T(6)«Sl*T(8)*CI)/C-WUL*T(2)*Sl)/T0(J) 

7  F(3)«(  OM#ULM*(T(B) *(APP— I*0)*C#T (2) *S1 ) ♦ (T ( 1 ) ♦ (T (5)*C1- 

i  ,♦  T(7)*Sl)/C)#DFG(J) *2*0#FG( J>* 

♦  ( T ( 5 ) *S1 ♦ T (7)#Cl)/C*WUL*T(l) *CI) 

F  (4)«(-0M*ULM*(T (7) ♦ (APP-l »0) *C*T (1)*S1) * (T (2) *(T (6)*Ct- 

♦  T(8)»S1)/C)*OFG(J)*2.0»FG(J)* 

♦  (T(6)*S1*T (8)*CT)/C*WUL*T(2)*C1) 

F  (5)  «-OM*T  (2)  *T  (7) « (PHO  ( J*  1 )  -PHO  ( J) )  /DS 
F(6)«  OM*T(l)*T(8)*(PHO(J*l)-PHOU))/DS 
F (7)«-T (5)* (PHO ( J*1)-PH0 ( J) ) /DS 
F  (8)  »>  T  (6)  *  (PHO  ( J*1 ) -PHO  ( J) ) /DS 
RETURN 

50  F(1)>0.0 

F (2)«0.0 
i  GO  TO  7 

END 
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SURFACE  SHIP  IN  A  SEAWAY  TOWING  A  DEEPLY  SUBMERGED  WEIGHT 
OMWAY 
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STEADY-STATE  SOLUTION 
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•i 

I 

1 

i 


CXM(J) 

DCXU) 

CYH(J) 

OCY(J) 

J 

0.000000 

0.000000 

0.000000 

0.000000 

1 

.007537 

-2.303692 

•012449 

.637901 

2 

.056209 

-1.564765 

•014046 

1.576627 

3 

•Q49291 

-.801759 

•017611 

2.339634 

4 

.051531 

-.302560 

.023934 

2.839033 

5 

•052516 

-.034669 

•030007 

3.106924 

6 

•050056 

.061509 

.033937 

-3.0600H4 

7 

.064777 

.079939 

.035120 

-3.061654 

6 

•037996 

-.034296 

•033633 

3.107297 

9 

.031219 

-.276547 

.031096 

2.865046 

10 

•026169 

-.671205 

•028626 

2.470387 

11 

•024713 

-1.195132 

.029620 

1.946461 

12 

•027646 

-1.725520 

.036254 

1.416073 

13 

•033999 

-2.166766 

.046134 

.972627 

14 

•042539 

-2.530422 

•064647 

.611171 

15 

.052629 

-2.841067 

.085431 

.300525 

16 

•064073 

-3.121891 

.110616 

.019702 

17 

•076691 

2.898573 

.140647 

-.243020 

16 

•091202 

2.647597 

.176153 

-.493996 

19 

•107182 

2.404798 

.217931 

-.736794 

20 
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APPENDIX  C 

LISTING  AND  SAMPLE  OUTPUT  OF  OMFLO 


PROGRAM  OMFLO ( INPUT (OUTPUT « TAPES- INPUT t TAPE6-0UTPUT ) 

DIMENSION  PHO  (200)  t  TO (200) (  U0(402)(  VO (200)  (  PHRA(200)t 

♦  PHRB(200) (  PHIA(200) *  PNlB(200)t  TRA(200>(  TRB(200).  TIA(200)« 

♦  TIB(200) (  URA (200) t  URB(200)«  UlA(200).  UIB(200)*  VRA(200)t 

S  ♦  VRBI200) t  VIA (200) t  V1B(200) (PMR(200) t  PHI (200).  TR(200>( 

♦  TI(200).  UR (200) t  UI(200).  VR(200).  VI (200).  FL(200).  FG(200)( 

♦  OFL(200) (  DFG(200) .  T(8)t  QPMO(200)t  QT0(200)t  OPHRA(200)t 

♦  QPHRB(200) f  OPHIA (200) .  OPHIB(200).  OTRA(200).  QTR8(200). 

♦  OTIAI200) (  QTIB(200) .  QURA(200)»  QURB(200)t  QUIA(200). 

10  ♦  QU1B(200) .  QVRA(20C) »  QVRB(200)t  QVIAI200).  QVIB(200) 

COMMON  DS(WUL(D(ALAM0(ALAM1 .ALAM2.BLAM1 . 

1BLAM2. AGAMO . AGAM 1 . AGAM2 . BGAM 1 » BGAM2  * 
20M(ULM(APP(C(PH0(DFL(DFG*K. I (TO.PHRtPHI iFL.FG 
EXTERNAL  -OAUXtOARN 
15  RE  AO  (Sdl)  T1TLE( 

♦  BL0(  BD0(  BV( 

♦  BLM(  DBL(  BDM(  DBD.  OM, 

♦  BM(  BMS.  BMH( 

♦  WUL(  ULM(  OS(  N( 

20  ♦  C(  ON(  CH(  CD(  AMP( 

♦  ALAMO (  ALAM1(  ALAM2(  BLAM1 (  BLAM2( 

♦  AGAMO (  AGAMl (  AGAM2.  BGAM1(  6GAM2 
C  WRITE  TITLE  ANO  INPUT  DATA 

WRITE  (6(21)  TITLE 

25  WRITE  (6(22)  BL0(  BOO.  BW 

WRITE  (6(23)  BLM(  DBL.  8DM(  DBD.  OM 
WRITE  (6(24)  BM(  BHS.  8MH 
WRITE  (6(25)  WUL(  ULM(  OS.  N 
WRITE  (6(26)  C(  ON(  CH(  CO(  AMP 

30  WRITE  (6(27)  ALAMO (  ALAM1(  ALAM2(  BLAM1(  BLAM2 

WRITE  (6(28)  AGAMO (  AGAMl (  AGAM2(  BGAMl.  BGAM2 
C  INITIALIZE 

0  -  0.5*DN*C*C*CM*CD 
AMC  ■  3.»4159*CH*CH«»DN»ANP/4.0 
35  APP  •  1.0  ♦  AMC/ULM 

CXO  -  0.0 
CYO  ■  0.0 
CXR  ■  0.0 
CXI  -  0.0 
40  CYR  -0.0 

CYI  ■  0.0 
KA  ■  N-l 
IA  «  N-l 

C  SOLVE  STEADY-STATE  PROBLEM 
45  ERR-.00001 

ERA*. 001 

CALL  ATA(BOO(YK) 

PHOB  -  ATAN2(BL0-BW(  -YK) 

TOB  ■  SORT (BD0*BD0  ♦(BL0-BW)*(BL0-BW) ) 

50  WRITE  (6(31) 

WRITE  (6(32) 

PHO ( 1 )  -  PHOB 
T0(1>  ■  TOB 

C  MULTIPLY  TO  FINO  VELOCITY  COMPONENTS 
55  DO  101  K  ■  1(N 


35 


UO  (K)  >  C*SIN(PHQ(K) ) 

VO ( K >  ■  C*COS(PHO (K) ) 

J  ■  -  K 

WRITE  (6*33)  PHO(K) »  TOCK)*  UO(K>*  VO(K>*  CXO*  CYO*  J 
60  IF  (K.EQ.N)  60  TO  101 

C  INTEGRATE  TO  FINO  TENSION*  ANGLE*  ANO  CABLE  CONFIGURATION 
T(1)-T0(K) 

T(2>-PHOOO 
TO)  ■  CXO 
65  T (A)  •  CYO 

FIRST*. 0 
ALS»FLOAT (K)*OS 

CALL  KUTMEH (A* ALS»T*ERR*OS*FIRST *HCX*ERA*OAUX) 

T0(K*1)>T(1) 

70  PH0(K*1)«T(2) 

CXO  *  T (3) 

CTO  *  T (4) 

101  CONTINUE 

C  SOLVE  DYNAMIC  PROBLEM 
75  WRITE  (6*34) 

WRITE  (6*35) 

C  REVERSE  LISTS  OF  PHO  ANO  TO 
DO  301  1  ■  1*  N 

QPHO(I)  ■  PHO(I) 

60  QTO ( I )  *  TO  ( I ) 

301  CONTINUE 

00  302  I  «  1*  N 

L  «  N*l-I 

PHO (I)  *  OPHO(L) 

85  TO ( I )  *  OTO(L) 

302  CONTINUE 

C  COMPUTE  VISCOUS-FORCE  TERMS  ALONG  CABLE 
DO  201  I  *  1*  IA 

FL ( I )  *  0* (*ALAMO“ALAMl*COS(PHO(I) ) ♦ALAM2*COS(2.0*PH0 (I) ) ♦ 
90  *  9LAM1*SIN (PHO ( I ) )-BLAM2*SIN (2.U*PH0 (I) ) ) 

FG( I )  ■  0* ( -AGAMO ♦ AGAM 1 *COS (PHO ( I ) ) -AGAM2*C0S ( 2 . 0*PH0 ( I ) ) 

♦  -8GAM1»SIN (PHO (I) ) ♦BGAM2*SIN (2.0*PHO (I) ) ) 

OFL(I)  *  D*(*ALAM1»SIN(PHO(I))-2.0*ALAM2*SIN(2.0*PHO(I)) 

♦  ♦BLAM1*COS(PHO(I) )-2.0»BLAM2*COS(2.0*PHO(I) ) ) 

95  OFG(l)  ■  0* ("AGAMl*SIN (PHO ( I ) ) ♦2.0*AGAM2*SIN(2.0*PH0 ( I ) ) 

♦  -BGAM1*C0S (PHO (I) ) *2. 0*8GAM2*C0S (2. 0*PH0 (I) ) ) 

201  CONTINUE 

C  COMPUTE  MODAL  SOLUTIONS 
C  COMPUTE  A  -MOOE  SOLUTION 

100  C  SET  VELOCITY  COMPONENTS  EQUAL  TO  ZERO  AT  BOTTOM 

URA(l)  *  0.0 

UIA(l)  *  0.0 
VRA(l)  >  0.0 
VI A ( 1 )  b  0.0 

105  C  ASSIGN  VALUES  TO  ANGLE  ANO  TENSION  AT  BOTTOM 

PHRA(l)  *  o.ooi 
PHI A ( 1 )  >  0.0 

TRA(l)  b  0.0 
T 1 A ( 1 >  b  0.0 

110  C  FINO  VELOCITY  COMPONENTS*  ANGLE*  ANO  TENSION  BY  INTEGRATING  UP 
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115 


120 


125 


130 


135 


140 


145 


150 


155 


160 


165 


C  THE  CABLE 

DO  202  1  ■  1*  1A 

T(l)-PHRAU) 

T (2)*PHIA(I> 

T(3)*TRA(I) 

T(4)-TIA<1> 

T(5)«URA<I> 

T (6)-UIA(I) 

T (7)»VRA(I) 

T(8)-VIA(I> 

FIRST-. 0 
ALS-FLOAT ( I ) *DS 

CALL  KUTMER (8* ALStT tERR*OS*FIRST t HCX*ERA»DARN) 

PHRA<IM)-T<1> 

PHIA<I*1>-1 (2) 

TRA(l*l)«T (3) 

TIA(I«1)-T(4) 

URA(I*1)-T (5) 

U1AU«1)«T<6> 

VRA  C 1*1 ) -T (7) 

VIA(I«1)-T<8) 

202  CONTINUE 

C  COMPUTE  B  -MODE  SOLUTION 

C  SET  VELOCITY  COMPONENTS  EOUAL  TO  ZERO  AT  BOTTOM 
URB(l)  «  0.0 
UIBd)  «  0.0 
VRBtl)  -  0.0 

VIB(l)  «  0.0 

C  ASSIGN  VALUES  TO  ANGLE  ANO  TENSION  AT  BOTTOM 
PHRB(l)  »  0.0 

PHIBUI  ■  0.0 

TRB(l)  -  0.1 
TIB ( 1 )  -  0.0 

C  FIND  VELOCITY  COMPONENTS*  ANGLE*  ANO  TENSION  BY  INTEGRATING  UP 

C  THE  CABLE 

DO  203  I  ■  1*  IA 
T ( 1 > -PHRB ( I ) 

T <2)-PHIB(I> 

T (3)-TRB(II 
T (4) -TIB (I) 

T (S)-URB(I) 

T(6)-UIB<I) 

T(7>-VRB(I) 

T(B)-VIBU) 

FIRST*. 0 
ALS-FLOAT U)*DS 

CALL  KUTMER (8* ALS*T»ERR*DS*FIRST»MCX*ERA*DARN) 

PHRR(IM)-TU) 

PHIB(I*1 )-T(2) 

TRB«I*1»-T(31 
TIB ( 1^1 >-T  <4> 

URB(IM)-T<5> 

UIB<IM)-T(6> 

VRBU«1)>T(7) 

VIB(Itl»-T(8) 
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203  CONTINUE 
C  REVERSE  LISTS 


OF  PHO  •  TO  *  A  NO  MODAL  SOLUTIONS 


00  303 

K 

■  It  N 

QPHO(K) 

■ 

PHO (K) 

170 

QTO(K) 

m 

TO(K) 

OURA(K) 

m 

URA(K) 

OUIA(K) 

u 

UIA(K) 

OVRA(K) 

* 

VRA(K) 

QVIA(K) 

m 

VIA(K) 

175 

QPHRA(K) 

•  PHRA(K) 

QPHIA(K) 

«  PH1A(K) 

OTRA(K) 

m 

TRA(K) 

0T1A(K) 

m 

TIA(K) 

QURB(K) 

m 

URB(K) 

180 

QUIB(K) 

9 

UIB(K) 

QVRB(K) 

9 

VRB(K) 

OVIB(K) 

9 

VIB(K) 

QPHRB(K) 

t 

■  PHRB(K) 

OPHIB(K) 

t 

■  PHIB(K) 

185 

QTRB(K) 

9 

TRB(K) 

QTIB(K) 

9 

TIB(K) 

303 

CONTINUE 

DO  304 

K 

■  1*  N 

M  ■  NM-I 

< 

190 

PHO (K1 

■ 

OPHO(M) 

TO (K)  ■ 

i  QTO(M) 

URA(K) 

9 

QURA(M) 

UIA(K) 

9 

QUIA (M) 

VRA(K) 

9 

QVRA (Ml 

195 

V1A(K) 

9 

QV1A(M) 

PHRA(K) 

9 

QPHHA(M) 

PHIA(K) 

9 

QPHIA(M) 

TRA(K) 

9 

QTRA(M) 

TIAdO 

9 

QTIA(M) 

200 

URB(K) 

9 

QURB (Ml 

UIB(K) 

9 

OUIB(M) 

VRB(K) 

9 

QVRB(M) 

VIB(K) 

9 

OVIB (Ml 

PHRB(K) 

9 

. OPHRB(M) 

205 

PH1B(K) 

9 

QPHIB(M) 

TRB(K) 

9 

QTRB(M) 

TIB(K) 

9 

QTIB(M) 

304 

CONTINUE 

210 


215 


220 


C  SUPERIMPOSE  MODAL  SOLUTIONS 

C  COMPUTE  FORCE  COMPONENTS  AT  TOP 


BOR 

BOI 

BLR 

BLI 

PR 

PI 

OR 

01 

PAR 


BOM*COS(BOB) 

BDM*5!N(BDB) 

BLM*COS(DBL> 

BLM*SIN(DBL) 

-  BDR*COS (PHO (111 

-  BOI#COS(PHO(1) ) 
BOR#SIN(PHO(1) )  ♦ 
BDI*SIN(PH0(1)>  * 


♦  BLR*SIN (PHO  Cl)> 

♦  BLI*SIN(PH0(11 1 
BLR*COS (PHO ( 1 ) ) 
BLI*COS(PHO(1) ) 


•  T*»*fl)  -  OM#(*(BMH«BMS1*UIA(11#COS(PHO(111#SIN(PHO(11 1 


♦  -(BM*BM  OS (PMO (11) *COS (PHO (11) ♦BMH#SIN (PHO (1)1 #SIN (PHO (!)))• 
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-  ■»*«*- uttfc 


2*5 

230 

235 

240 

245 

250 

255 

260 

265 

270 

275 


♦  man 

PA!  ■  TIA(l>  -  OM«t-(BMH-BMS)4URAU)*COStPHOUn«SIN(PHOCin 

♦  ♦ (BM*BMS*COS  (PHO ( 1 ) ) *COS (PHO (1 ) ) ♦BNH*SIN (PHO ( 1 > ) «SIN (PHO (1) ) ) • 

♦  VRAUM 

PSR  ■  TRB(i)  -  OM*K(BMH-BMS>*UlB(I)*COS(PMO(in«SIN(PHO(in 

♦  - (BM*BMS*C0S (PHO (1 ) ) *COS (PHO ( 1 ) ) ♦BHH*SIN (PHO ( I )) *SIN (PHO (1 )>  > * 

♦  VIBID) 

PBI  ■  7XBT1)  •  OM»(-(BMH-BMS)*URB(l>*COS(PHO(l>)*SIN(PHO(l)> 

♦  ♦  (BM*BMS*COS  (PHO  a  n  *COS  (PHO  ( n  )  ♦BMH*SW  (PHO  ( 1 ) )  •SIN  (PHO  ( 1 )  > )  • 

♦  VRB(1>> 

OAR  ■  TO(I)*PHRA(I)“ON*(*  (BM«BMS*SIN(PH0(1))#SIN(PH0 (1)  I  ♦ 

♦  6MH*COS(PHO(l>)*COS(PHO(l))>*UIA(l)-(BMH-eMS)*VIA(I)*COS(PHO(l>> 

♦  •SIN(PH0(1>>> 

OA!  ■  TO ( I)*PHIA( 1 )-0M*(-(BN*8MS*SIN (PHO ( 1 ) )#SIN(PHO ( 1 ) I  ♦ 

♦  BMH«COS<PHO(in*COS(PHO(in>*UHA(l>*(BMH-BMS)*VRA(l>*COS(PHO(in 

♦  •SIN(PHOUn) 

OBR  ■  TO<l)*PHRB(n-OM*(«(BH«BMS*SIN(PHO(in*SIN(PHO(in* 

♦  BMH*COS<PHO<1))»COS(PHO(I)>)»UIB(1)-(BMH-BMS)*VIB(1>*COS(PHO(1)) 

♦  •siN(pHoan) 

OBI  «  TO  ( 1 )  #PHIB  ( I )  “OH*  (- (BM*BMS*SIN  (PHO  ( 1 ))  *SIN  (PHO  ( 1  >  >  ♦ 

♦  BMH*COS  (PHO  (I))  *COS (PHO  ( 1 )  >)  •URB(  1 )  ♦ (BNH-BMS)  *VRB  ( I )  *COS  (PHO  ( I )) 

♦  •SIN(PH0(1)>) 

DENOMR  ■  P AR*QBR-P A I *081 -PBR*UAR*PB  I  *0A  I 
OENOHI  ■  PA  I  #QBR  *P AR*OB  I  -PB  I  *UAiR-PBR*OA  I 
ANUHR  ■  P  R*QBR-P  I«QBI-PBR*0  R*PBI*Q  I 

ANUMI  ■  P  I*OBR*P  R#0BI“PBI*O  R-PBR*Q  I 

BNUHR  •  PAR*Q  R-PAI*Q  I-P  R*UAR*P  1*0 A I 

BNUMI  ■  PAI*0  R*PAR*Q  I-P  I*UAR-P  R*QAI 

0EN0H2  «  OENOHR*OENOMR  ♦  OENOHI •OENOHI 
CAR  ■  ( ANUHR *OENO MR  ♦  ANUMI •OENOHI )/D£N0M2 

CAI  ■  (-A'(UHR •OENOHI  ♦  ANUH I *DENOMR )  /0EN0H2 
CBR  -  (BNUHR •DENOMR  ♦  BNUH I* OENOHI )  /DEN0H2 
CBI  ■  (-BNUHR  •OENOHI  ♦  BNUH  I  *DENOHR )  /DEN0H2 
C  COMPUTE  VELOCITY  COMPONENTS*  ANGLE*  AND  TENSION  ON  CABLE 
DO  204  K  ■  1*  N  . 

UR  (K)  ■  CAR*URA (K)  •  CAI»U1A(K)  ♦  CBR*URB(K)  -  CBI*UIB<K) 

UKK)  ■  CAI*URA(K)  ♦  CAR*UIAOO  *  CBI*URB(K)  ♦  CBR*UIB(K) 

VR(K>  «  CAR*VRA(K)  -  CAI*VIA(K)  *  CBR*VRB(K)  -  CBI*VIB(K) 

VI (Ki  ■  CAI*VRA (K)  *  CAR*VIA(K)  *  CBI*VRB(K)  ♦  CBI*VIB(K) 

PHR(K)  •  CAR*PMRA(K)  -  CAI«PHIA(K)  *  CBR*PHRB(K)  -  CBI*PH1B(K> 

PHI (K)  ■  CAI*PHRA (K >  ♦  CAR*PH1A(K)  ♦  CBI*PHRB(K)  ♦  CBR*PHIB(K) 

TR(K)  ■  CAR*TRA(K)  -  CAI»TIA(K)  ♦  CBR*TRB(K)  -  CBI*TIB(K) 

TI (K)  «.  CAI#TRA (K)  *  CAR*TIA(K)  *  CBI*TRB(K)  ♦  CBR*TIB«K) 

C  COMPUTE  PHYSICAL  VELOCITY  COMPONENTS  ON  CABLE 
URP  ■  C#PHR(K)*COS(PHO(K) )  ♦  UR(K> 

DIP  ■  C*PHI  (K)*COS(PHO(K) )  *  UKK) 

VRP  ■  -  C*PHR(K)*SIN(PHO(K) )  ♦  VR(K> 

VIP  ■  -  C*PHI (K) *SIN (PHO (K) )  *  VI(K> 

C  COMPUTE  MAGNITUDE  ANO  PHASE  OF  VELOCITY  COMPONENTS*  ANGLE*  A  NO 
C  TENSION  ON  CABLE 

UM  ■  SORT  (URP*URP  ♦  UIP*UIP) 

CALL  ATA(URP»  YK) 

OU  ■  ATAN2(UIP*  YK) 

VM  ■  SORT  (VRP* VRP  ♦  VIP*VIP) 

CALL  ATA(VRP»  YK) 
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280 


285 


290 


295 


300 


305 


310 


315 


320 


325 


330 


0V  *  ATAN2 (VIP*  YK) 

PHM  *  SUKT<PHH<K)*PHR(K>  ♦  PHI (K)*PHI (K) ) 

CALL  ATA(PHR(K)  t  YK) 

DPH  »  ATAN2(PH1 (K) *  YK) 

TM  a  SQKT <TR(K)*TH(K>  ♦  TI(K)*TI(K)) 

CALL  ATA( TR(K) »  YK) 

OT  >  ATAN2ITI < K ) *  YK) 

J  «  -  K 

WRITE  (6*36)  PHM*  DPH*  TH*  DT*  UM*  DU*  VM»  DV*  J 

204  CONTINUE 
WRITE  (6*37) 

C  COMPUTE  MAGNITUDE  ANO  PHASE  OF  DYNAMIC  CABLE  CONFIGURATION 
DO  2 05  K  =  1 *N 
CXM  *  SUWT<CXH*CXK*CXI«CXI> 

CALL  ATAiCXR*  YK) 

OCX  *  ATAN2 (CX I  *  YK) 

CVM  *  SORT  <CYR*CYR*CYI*CYI ) 

CALL  ATATCYR.  YK) 

DCV  *  ATAN2TCYI «  YK) 

J  «  -  K 

WRITE  (6*38)  CXM*  OCX*  CYM*  DCY*  J 
IF  (K.EQ.  N)  GO  TO  205 
CXR  a  PHR(K)«SIN(PHO(K))*DS 
CXI  *  PHI (K)*SIN(PH0(K) )#0S 
CYR  ■  -PHH  <K) #COS (PHO (K) ) *DS 
CYI  »  -PHI <K)*COS<PHO<K)>*D$ 

205  CONTINUE 

11  FORMAT  (  10X*  A6  /  10X*.  3F10.5/  10X*  5F10.5/  10X*  3F10.5/ 

♦  lOXt  3F10.5*  110/  10X*  5F10.5  /  10X.  5F10.5/  10X*  5F10.5) 

21  FORMAT  (1H1*  10X*  S3HDEEPLY  SUBMERGED  SUBMARINE  TOWING  A  FLOAT  IN 

♦A  SEAWAY  /  1H-*  13X*  A6///) 

22  FORMAT  (  IX.  8HBL0  •  .  F10.5*  4X*  8HBD0  a  »  F10.5<  4X* 

♦  BHBW  ■  .  FlO.b  //) 

23  FORMAT  (  IX*  8HBLM  *  •  F10.5*  4X*  8HDBL  a  *  F10.5*  4X. 

♦  SHBDM  >  *  F10. 5.  4X.  8H0BD  a  .  F10.5*  4X*  8H0M  ■  * 

♦  F10. 5  //) 

24  FORMAT  (  IX*  8HBM  a  *  F10.5*  4X*  8HBMS  a  ,  F10.5*  4X* 

♦  BHBMH  a  ,  F10.S  //) 

25  FORMAT  (  IX*  8HWUL  ■  *  F10.5*  4X*  8HULM  a  ,  F10.5*  4X* 

«  BHDS  a  ,  F10. 5*  4X*  8HN  a  ,  HO  //) 

S'S  FORMAT  (  IX.  8HC  a  ,  F10.S*  4X*  8HDN  a  *  F10.5*  4X* 

♦  8HCH  a  ,  F10. 5*  4X»  8HCD  '  »  »  F10.5*  4X.  8HAMP  a  .F10.5//I 

27  FORMAT  <  IX*  8HALAM0  a  ,F10.S*  4X.  8HALAM1  a  ,  F10.5.4X* 

♦  BHALAM2  a  ,  F10.5*  4X*  8HBLAM1  a  ,  F10.5*  4X*  8HBLAM2  a  « 

♦  F10.5  //) 

28  FORMAT  (  IX*  BHAGAMO  *  «  FlO.b*  4X*  SHAG AMI  a  «  F10.5*  4X* 

♦  8HAGAM2  a  ,  HO. 5*  4X*  8HBGAM1  a  ,  F10.S*  4X»  8HBGAM2  a  , 

♦  F10.5  //) 

31  FORMAT  (  22H1 STEADY-STATE  SOLUTION) 

32  FORMAT  (  1H0*  IX*  6HPH0 ( J) *  10X*  6HT0TJ)  *  10X*  6HU0TJ)  *  JOX* 

♦  6HV0 ( J)  «  10X*  6HCX0 ( J) *  10X*  6HCY0 ( J) *  10X*  1HJ  //) 

33  FORMAT  <  2X»  F9.6*  7X*  F13.6*  3X*  2(  F11.6*  5X).  2IF13.6*  3X) *  14) 

34  FORMAT  (  17H1DYNAMIC  SOLUTION) 

35  FORMAT  (  1H0.  IX*  6HPHM(J)»  9X»  6HDPHT J) »  9X*  6HTHTJ)  .  9X* 

♦  6HDTTJ)  •  9X»  6HUMTJ)  •  9X*  6HUUTJ)  *  9X*  6HVMI J)  «  9X.  6HDVT J)  * 
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♦  9X»  1HJ  //) 

36  FORMAT  <  2X»  2(  F9.6*  6X)»  F12.6*  3X»  F9.6*  6X.  2(F10.6«  5X. 

♦  F9.6*  6X1*  14) 

37  FORMAT  (  1H1  *  IX*  6HCXMU).  10X*  6H0CX ( J)  »  10X*  6HCYM(J>»  10X* 

335  ♦  6MDCY ( J) •  10X*  1HJ  //» 

36  FORMAT  (  2X.  2(  F12.6*  4X*  F9.6*  7X).  14) 

STOP 

END 


S 
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SUBROUTINE  ATA(X.YK) 

C  THIS  SUBROUTINE  PREVENTS  ATAN2  ¥  ROM  DIVIDING  BY  ZERO  MHEN 
C  COMPUTING  THE  ARCTANGENTS  OF  PI/2  AND  3*PI/2 
ER  ■  0.00001 

IF  (A8S(X) .GT.EN)  GO  TO  10 
YK  ■  ER 
RETURN 
10  YK  >  X 
RETURN 
ENO 
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SUBROUTINE  KUTMCR(N*T*VO*EPS*H*FIRST*HCX*A*OAUX> 

C  THIS  SUBROUTINE  PERFORMS  FOURTH-ORDER  KUTTA-MERSON  INTEGRATIONS 
DIMENSION  YO(10)*Y1(10)*V2(10)*FO(10)*F1MO)«F2(10) 

IF (FIRST) 20 *10*20 
10  HC*H 
IPLOC-1 
FIRST-1. 

20  IOC-0 

HCX-MC 

30  CALL  OAUX(T*YO«FO) 

39  DO  40  1*1  *N 

40  Y1(I>-Y0(I)«(HC/3.)*F0(I) 

CALL  0AUXU*HC/3.*Y1*F1) 

00  50  1*1. N 

SO  Y1(I>-Y0<I)*<HC/6.)*F0(I)«(HC/6.)*F1(I> 

CALL  DAUX<T*HC/3.«Y1*F1) 

00  60  I-1*N 

60  V1(I)-Y0(I) *HC/8.*F0 (I) ♦»375*HC#F1 ( I ) 

CALL  0AUX(T«HC/2.«Y1*F2) 

00  70  1*1 *N 

70  YI  (I )*Y0( I 1  *HC/2.*F0 (I)-1.5*HC*F1 (I )  ♦2.#MC*F2(I) 

CALL  0AUX(T*HC*V1»F1) 

DO  60  1*1 *N 

00  Y2(I )*YO(I) *HC/6.*F0 (I) *.66666667*HC*F2( I) *(HC/6.)*F1 (I) 

INC-0 

00  110  1*1 *N 
ZZZ-ABS(Y1(I))-A 
IF(ZZZ)  8S.87.S7 
85  ERROR*ABS ( «2*( VI ( I ) "Y2 ( I ) ) ) 

IF(ERROH-A) 100*100*90 
87  ERROR*ABS( «2-.2*Y2 III/YKIII 
IF (ERROR-EPS) 100*100*90 

90  X*128.*ABS(HC)-ABS(H> 

IF (X)  91*95*95 

91  WRITE(6*92)  T. ERROR 

92  FORMAT (21M  RELATIVE  ERROR  AT  T-  1P1E12.3.3HIS  F10.6  > 

FIRST  ■  2. 

RETURN 
95  HC-HC/2. 

IPLOC-2  »IPLOC 
LOC-2  *LOC 
HCX-HC 
GO  TO  30 

100  IF(ERR0R*64.-EPS) 110*110*101 

101  INC*1 

110  CONTINUE 

111  T«T*MC 

00  112  1*1 *N 

112  Y0(I)*Y2(I) 

LOC-LOC+1 

IF (LOC-IPLOC) 120*210*210 
120  IF (INC >210* 130 *210. 

130  IF (L0C-(L0C/2) *2)210*140*210 
140  IF (IPLOC-1) 210*210*200 
200  HC*2.*HC 


KUTMOOSO 
KUTM0060 
KUTM0070 
KUTMOOSO 
KUTM0090 
KUTM0100 
KUTM01 10 
KUTM0120 
KUTM0130 
KUTM0140 
KUTM01SO 
KUTM0160 
KUTM0170 
KUTM01B0 
KUTM0190 
KUTM0200 
KUTM0210 
KUTM0220 
KUTM0230 
KUTM0240 
KUTM0250 
KUTM0260 
KUTM0270 
KUTM0280 
KUTM0290 
KUTM0300 
KUTM0310 
KUTM0320 
KUTM0330 
KUTM0340 
KUTM0350 
KUTM0360 
KUTM0370 
KUTM0390 
KUTM0390 
KUTM0400 
KUTM0410 
KUTM0420 
KUTM0430 
KUTM0440 
KUTM0450 
KUTM0460 
KUTM0470 
KUTM04B0 
KUTM0490 
KUTM0500 
KUTM0510 
KUTM0520 
KUTM0530 
KUTM0S40 
KUTM0550 
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1 


60 


5 


10 


15 


20 


25 


LOOLOC/2 

kUTMOSAO 

IPL0C-IPL0C/2 

KUTMOSTO 

210 

IF ( I PL OC-LOC) 30*220*30 

KUTMO50O 

220 

RETURN 

KUTM0S90 

END  , 

i 

i  i  ‘  | 

KUTM0600 

I 

, 

1  '  i  < 

SUBROUTINE  DAUX(Z*TX»F)  , 

THIS  SUBROUTINE  PROVIDES  THE  INTEGRANDS  FOR  KUTMER 

/ 

1 

WHEN  IT 

C  COMPUTES  THE  STEADY-STATE  TENSION*  ANGLE  *  AND  CABLE  CONFIGURATION 
DIMENSION  TX (8) «  F (B)«PHO (200) .DFL (200) ♦DFG<200) .TO (200) . 

♦  PHR (200) ♦  PHI (200)  1  , 

,  COMMON  DS*WUL*D*ALAM0*ALAM1*ALAM2*BLAM1* 
1BLAM2*AGAM0»AGAM1*AGAM2,BGAM1*BGAM2»  i  , 

20M*ULM» APP.C*PHO*DFL*DFG*K* I »  TO  tPHR.PHI 
TX1«TX(1) 

TX2«TX(2)  '  ,  , 

C1*C0S(TX2)  i 

S1*SIN(TX?)  , 

C2*C0S(2.*1X2)  ,  ,  1 

S2»S1N(2.*1X2)  ,  / 

FL*D# (ALAM0-ALAM1*C1*ALAM2*C2*BLAM1*S1-BLAM2*S2)  „ 

FG«D* (-AGAM0*AGAM1*C1-AGAM2*C2-«GAM1«S1*BGAM2»S2) 

IF(TX1 .EO.O.)  GOTO  50 

F(2)»-(WUL*C1-FL)/TX1  i  i 

F  (3)  =  -  Cl 

,  F  (4)  *  -  SI 

7  ,  FUJa-rfWUL-SUFG)  , 

RETURN 

50  F  (2)=0.0 , 

'  GO  TO  7  1 
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25 


30 
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i  SUBROUTINE  DARN(Z*T«F> 

C  >  THIS  SUBROUTINE, PROVIDES  THE  INTEGRANDS  FOR  KUTMER  WHEN  IT 
C  COMPUTES  THE  DYNAMIC  ANGLE*  TENSION*  AND  VELOCITY  COMPONENTS 
DIMENSION  T (81  *  F (8) *PH0 (200> tDFL(200) *DFG(200) *T0(200) * 

♦  FL1200) *FG<200> *PHR(200) »PHI (200) , 

COMMON  DS*MUL*D*ALAM0*ALAM1 *ALAM2»BLAM1 * 

.  1BLAM2*AGAM0»AGAM1*AGAM2*BGAM1»BGAM2* 

2OM»ULM*APP*C*PHQ«0FL»DFG*K»I*T0*PHR*PHI*FL*FG 

C1«COS(PHO(I)> 

Sl«SIN(PWni)) 

IF(TOd).EQ.O.O)  GO  TO  50 

F(l)»(-OM*ULM#(Ar*P*T (p) ♦  (APP-1.01*C*T (2)*C1)-T I3)*(PH0(I*1)-PH0 (I) 
♦)/D$-(T(l)*(f (5)*Cl-T(7)*Sl)/t)*DFL(I)  , 

♦  -2«0#FL (I)*(T(5)*SltT (7)*C1)/C-WUL#T(1)#S1)/T0(I) 
F(2)»<  OM*ULMP(APP*T(5)*(APP-1.0)*C*T(l),*Cl)-TU)*(PH0(I*n-PH0m 

♦  )/DS-(T(2)MT(6)*Cl-T(8)*Sl>/C) *DFL(I) 

♦  -2.0»FLU)*(T(b)*Sl»T(8>«Cl)/C-BUL*T(2)»Sl)/T0(I) 

7  F(3)»(  0M*ULM*(T(8) ♦<APP-1«0)*C*T I2)*SD  ♦  (T  (1)  MT  (5l*Cl- 

♦  T(7»»Sl)/C)«DFG(I)*2.0»FGtn*  < 

♦  (T(5)*S1*T (7)*C1J/C*BUL#T (1)#C1) 

F  (4) ■ (-OM*ULM* ( T  (7)  ♦  (APP-1  »0)  #C*T  ( 1 ) *S1 )  ♦  (T i(2)  ♦  (T  (6)  #C1- 

♦  T (8)*S1)/C)*DFG(I) ♦2»0*FG(I)*- 

♦  (T(6)*S1*T  C8)  *Cl)/C*BUL#T C2)*CO 

F  (5)*-OMPT  (2)  *T  (7)#(PH0  (I*l)-PHO  (I  >’)/DS 
F (6) ■  0M*T (11 ♦T  (8)*(PH0(I*1)-PH0(I) )/DS  , 

F(7)»-T(5)*(PHO(I*1)fPHO(I) )/DS 

F  (8)*»-T  (6> •  IPH0  <  1*1  >^PH0 (I >  1  . 

RETURN  .  . 

50  F(1)«0.0 

F(2)»0.0 

GO  TO  7  ,  , 

END 
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STEADY-STATE  SOLUTION 
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CXM(J) 

0. 


DCX(J) 

CYNIJ) 

DCY(J) 

J 

000000 

0.000000 

0.000000 

0.000000 

-1 

002306 

1.695351 

.000924 

1.695361 

-2 

000940 

2.359663 

.000656 

2.359663 

-3 

000421 

3.063075 

.000415 

3.063075 

-4 

000224 

-2.570909 

.000267 

-2.570909 

-5 

00013S 

-2.065947 

.000213 

-2.0659*7 

-6 

000005 

-1.693579 

.000161 

-1.6935 

-7 

000055 

-1.357275 

.000119 

-1.406266 

-6 

000035 

-1.163365 

.000065 

-1.163365 

-9 

000024 

-.939157 

.000066 

-.939157 

-10 

000017 

-.729205 

.000052 

-.729205 

-11 

000013 

-.557606 

.000045 

-.557606 

-12 

000011 

-.457430 

.000042 

-.461557 

-13 

000010 

-.415907 

.000042 

-.457120 

-14 

000010 

-.456233 

.000043 

-.534646 

-15 

000009 

-.534252 

.000046 

-.673476 

-16 

000010 

-.624314 

.000050 

-.653303 

-17 

000010 

-.707*60 

.000055 

-1.055456 

-16 

000010 

-.773746 

.000061 

-1.275290 

-19 

000011 

-.519057 

.000066 

•1.424606 

-20 
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